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Abstract 

This  research  delves  into  the  area  of  image  enhancement  for  the  visually  impaired. 
Binocular  macular  degeneration,  a  visual  impairment,  affects  many  Americans;  since  this 
condition  could  not  be  corrected  with  conventional  glasses,  the  literature  suggested  using  an 
enhancement  system  which  used  a  pre-emphasis  algorithm  to  enhance  the  input  image  for 
output  to  the  observer.  The  work  of  Dr.  Eli  Peli,  a  pioneer  in  the  field  of  image  enhancement, 
is  examined  and  reproduced.  Since  his  work  concentrated  mainly  on  frequency  analysis  of 
images,  the  bulk  of  this  research  involves  using  discrete  wavelet  analysis  to  augment  that  work. 
A  biorthogonal  wavelet  set  is  used  to  enhance  images.  That  same  wavelet  set  is  used  to  provide 
a  simulation  of  a  person’s  perception  of  an  image  based  upon  that  person’s  contrast  sensitivity 
function.  The  wavelet  enhancement  follows  Dr.  Peli’s  method  of  enhancement  using  a 
modified  pre-emphasis  filter.  The  wavelet  enhancement  results  are  similar  to  that  technique, 
with  the  most  notable  difference  being  that  wavelet  enhancements  tend  to  accentuate  the 
horizontal  and  vertical  details  of  the  image  more  than  the  spatial  frequency  concentric  filter 
techniques  outlined  by  Peli  and  others.  The  perception  simulations  use  wavelet  analysis  and 
are  based  on  the  work  accomplished  by  Peli  concerning  “local  band-limited  contrast,’’  with 
the  frequency  analysis  of  that  work  giving  way  to  a  wavelet  analysis;  the  results  here  are 
significantly  different  than  those  results  seen  in  the  Peli  simulation.  A  method  to  enhance 
images  based  on  suprathreshold  characteristics  of  the  observer  is  outlined.  Methods  in  the  past 
have  concentrated  on  contrast  threshold  analysis,  so  the  suprathreshold  method  of  enhancement 
is  the  next  logical  step.  A  method  for  optimizing  localized  image  enhancement  is  presented. 
This  method  has  been  used  to  increase  the  resolution  for  pilots  viewing  virtual  environments 
and  could  take  advantage  of  the  existing  retinal  tracking  hardware  and  software  already  in 
existence  to  better  enhance  images  for  the  visually  impaired. 
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Image  Perception  Wavelet  Simulation  and  Enhancement  for  the  Visually 

Impaired 


/.  Introduction 

This  research  focused  on  enhancing  images  for  persons  suffering  from  binocular  macular 
degeneration,  or  loss  of  vision  in  the  foveal  area  of  both  eyes.  The  macula  lutea  contains  the 
highest  concentration  of  cones  (nerve  endings  giving  rise  to  color  vision)  on  the  retina;  the 
fovea  lies  at  the  center  of  the  macula  lutea  where  the  lens  focuses  incoming  light,  and  a 
person’s  visual  acuity  is  directly  correlated  to  the  condition  of  the  cells  in  this  area.  Persons 
suffering  from  central  scotoma  (darkness  in  the  center  of  the  eye)  due  to  cell  damage  can  have 
difficulty  reading  text  and  recognizing  faces.  Different  degrees  of  macular  damage  require 
different  visual  enhancement  techniques  so  people  with  low  vision  due  to  central  scotoma  can 
regain  some  of  their  original  visual  ability. 

Conventional  glasses  are  most  useful  if  a  person’s  visual  impairment  is  due  to  aberrations 
of  the  lens  of  the.  eye;  therefore,  if  people  suffering  from  central  scotoma  are  to  be  helped,  a 
different  method  of  correction  must  be  made  available.  Johns-Hopkins  is  currently  testing  a 
new  device  known  as  the  Low  Vision  Enhancement  System  (LVES,  pronounced  Elvis)  which 
consists  of  a  pair  of  binocular  cameras,  a  single  magnifying  camera,  and  two  overlapping  LCD 
screens  to  present  a  virtual,  enhanced  world  to  a  low-vision  patient  in  the  hopes  that  computer 
digital  enhancement  may  improve  the  patient’s  ability  to  recognize  faces  and  read  text. 

The  largest  part  of  this  thesis  effort  concentrated  on  helping  the  visually  impaired  with 
two  new  approaches,  both  using  a  biorthogonal  wavelet  set.  First,  the  wavelet  decomposition 
detail  coefficients  were  manipulated  based  on  a  viewer’s  contrast  sensitivity  function  in  order 
to  gain  insight  into  what  is  actually  perceived  by  a  patient  when  viewing  an  image.  The  choice 
of  wavelets  as  a  human  perception  diagnostic  tool  followed  from  the  work  of  Wang  (35)  in 
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which  he  modeled  human  auditory  perception  using  a  wavelet  basis  set.  The  premise  that 
human  corticle  cells  exhibit  some  of  the  characteristics  of  a  wavelet  processor  is  supported 
by  the  work  of  Jones  and  Palmer  (13)  concerning  modelling  of  simple  receptive  fields  in  the 
cortex  of  cats.  Second,  the  effects  of  boosting  certain  detail  coefficients  in  a  given  image  were 
explored  in  the  hope  that  this  manipulation  could  be  beneficial  to  a  low-vision  viewer. 

A  new  method  using  suprathreshold  characteristics  of  a  patient  for  enhancement  was 
also  outlined  as  was  using  only  the  information  in  an  image  in  the  area  of  interest  to  the 
observer  for  enhancement. 

1.1  Background 

Image  enhancement  for  the  visually  impaired  has  been  studied  since  1984  when  Dr. 
Eli  Peli  (26)  conjectured  that  digital  techniques  could  be  used  to  improve  face  recognition 
ability.  His  subsequent  research  explored  ways  to  enhance  facial  images  using  digital  image 
processing,  and  he  used  frequency  analysis  as  his  primary  tool  for  image  manipulation. 
Wavelets  have  been  used  extensively  as  a  method  of  extracting  information  contained  in 
images  (4,  14,  17,  18,  29).  Wavelet  decompositions  are  analogous  to  (but  should  not  be 
confused  with)  frequency  filter  analysis  and  their  implementation  during  image  processing 
parallels  the  analysis  done  by  Peli  concerning  both  image  enhancement  and  human  perception 
of  images. 

1.2  Problem  Statement 

Investigate  frequency  and  wavelet  analysis  in  conjunction  with  tested  individual  per¬ 
ception  parameters  to  simulate  the  view  of  an  image  to  a  person  with  macular  degeneration 
and  to  enhance  images  to  be  presented  to  persons  with  low  vision  due  to  binocular  macular 
degeneration. 
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1.3  Research  Questions 


Does  wavelet  analysis  lend  itself  to  a  better  understanding  of  how  people  with  visual 
defects  perceive  images?  Can  discrete  wavelet  analysis  be  used  for  effective  image  enhance¬ 
ment? 

1.4  Standards 

Dr.  Eli  Peli  pioneered  the  idea  of  using  digital  image  processing  to  enhance  images  for 
viewers  with  visual  defects  (26).  He  conjectured  these  enhancements  would  allow  low-vision 
viewers  to  more  readily  discriminate  details  about  their  visual  environment,  and  that  in  turn 
would  improve  their  ability  to  recognize  faces.  He  also  provided  the  concept  of  “band-limited 
contrast”  which  described  the  contrast  in  a  given  image  locally  and  allowed  experimenters 
to  simulate  what  a  person  with  measured  visual  parameters  might  observe  when  looking  at 
an  image.  He  observed  the  effects  of  various  enhancement  techniques  (outlined  in  Chapter 
II)  and  initiated  testing  of  two  techniques  (adaptive  filtering  and  adaptive  thresholding)  on 
patients  with  low  vision.  He  and  others  attempted  to  tune  the  image  enhancement  for  a  specific 
patient  by  using  the  contrast  sensitivity  function  of  that  patient  as  the  guide  for  enhancement. 
Their  enhancements  were  applied  to  the  entire  image. 

1.5  Scope 

This  research  first  investigates  and  reproduces  the  simulation  of  image  perception  based 
on  Frequency  analysis  presented  by  Peli  (23,  19).  With  those  simulations  reproduced,  the 
discrete  wavelet  transform  is  used  in  an  analogous  manner  to  simulate  human  perception. 
Images  are  enhanced  using  techniques  established  by  Peli  (26,  22,  20,  24,  21,  7)  and  Law- 
ton  (15)  and  then  compared  to  images  enhanced  using  discrete  wavelet  analysis  to  parallel  the 
aforementioned  work  by  direct  manipulation  of  the  decomposition  detail  coefficients  rather 
than  the  manipulation  of  the  spatial  frequency  spectrum.  After  analyzing  this  information,  a 
method  of  enhancement  based  on  measured  suprathreshold  characteristics  of  an  observer  is 
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suggested.  Finally,  recommendations  are  made  on  the  use  of  a  retinal  tracking  device  now 
available  for  enhancement. 

1.6  Methodology 

Since  this  thesis  was  the  first  at  AFIT  to  be  concerned  with  image  enhancement  for 
the  visually  impaired,  the  first  step  was  to  discover  what  work  had  been  done  in  that  area. 
A  literature  search  revealed  the  bulk  of  the  research  belonged  to  Dr.  Eli  Peli  of  the  Eye 
Research  Institute.  First,  his  visual  degradation  models  simulating  central  scotoma  (19) 
were  reproduced.  Second,  with  these  human  perception  simulations  duplicated,  Peli’s  face 
enhancement  techniques  (26,  22,  20,  24,  21,  7)  and  Lawton’s  letter  enhancement  techniques 
(15)  (as  applied  to  faces)  were  reproduced.  Third,  using  Peli’s  work  as  a  foundation,  a  new 
model  was  developed  using  a  biorthogonal  wavelet  set  to  simulate  image  perception  of  a 
person  with  a  given  contrast  sensitivity  function.  Fourth,  that  same  biorthogonal  wavelet  set 
was  used  to  enhance  images  based  on  a  given  person’s  visual  perception.  Fifth,  a  method  was 
devised  to  enhance  an  image  using  information  based  on  not  only  measured  contrast  threshold 
characteristics  of  an  observer  (as  was  done  in  the  past),  but  also  on  measured  suprathreshold 
characteristics  of  the  patient  combined  with  the  rms  contrast  of  the  image  being  enhanced. 
Sixth,  a  scheme  for  enhancing  only  the  area  where  the  patient  is  looking  using  a  retinal  tracking 
device  is  outlined.  Wavelets  have  never  been  used  in  the  context  of  benefiting  persons  suffering 
from  low  vision.  The  bulk  of  this  research  explored  the  value  of  wavelets  as  enhancement 
tools,  and  also  the  validity  of  using  wavelet  analysis  to  simulate  a  person’s  vision  based  on 
that  person’s  measured  contrast  threshold  perception  characteristics. 

1.7  Overview  of  Thesis 

Chapter  2  provides  a  background  of  digital  image  enhancement  methods  and  results  to 
date,  as  well  as  a  review  of  practical  application  of  wavelet  analysis  of  images.  Chapter  3  covers 
various  methodologies  used  to  enhance  images  and  simulate  visual  perception.  Chapter  4 
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discusses  the  results  of  these  previous  experiments.  Chapter  5  presents  concluding  arguments 
and  work  which  could  be  done  to  further  the  research  in  this  area. 
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11.  Literature  Review 


2.1  Introductton 

In  1972,  Stockham  investigated  some  of  the  mechanics  of  how  visual  information  is 
processed  in  humans  (32).  In  that  paper,  he  explored  how  the  human  visual  system  (HVS) 
processes  visual  information  and  pointed  out  that  our  view  of  the  world  is  not  the  actual  world, 
but  instead  a  simulation  seemingly  designed  to  maximize  the  amount  of  useful  information 
about  the  scene  being  viewed.  Ginsburg’s  1978  dissertation  concentrated  on  explaining  the 
HVS  in  terms  of  filtering  characteristics  alone  (8).  There,  he  conjectured  that  the  HVS  could 
be  modelled  as  octave  bandwidth  spatial  filters  with  varying  orientations.  If,  in  fact,  the  HVS 
does  behave  as  a  filter  bank,  then  measurement  of  the  frequency  response  of  the  system  would 
be  the  natural  first  step  in  any  subsequent  analysis  of  human  perception. 

Jones  and  Palmer  demonstrated  that  invasive  types  of  measurements  are  feasible  when 
analyzing  cats  (12),  but  invasive  techniques  are  not  usually  possible  with  human  subjects. 
Contrast  sensitivity  measurements  describe  the  minimum  contrasts  necessary  for  a  tested 
subject  to  perceive  a  set  of  given  spatial  frequencies.  These  measurements  allow  us  (the 
experimenters)  to  non-invasively  glimpse  what  may  be  perceived  by  another  person.  Since  this 
Thesis  concerns  itself  with  image  enhancement  for  low-vision  patients  suffering  from  binocular 
macular  degeneration,  knowledge  about  how  these  people  perceive  their  environment  could 
be  used  to  make  their  world  more  recognizable.  Although  presenting  the  world  through  an 
electronic  headset  would  not  be  the  same  as  seeing  with  perfect  eyes,  Stockham  already  pointed 
out  that  natural  human  visual  perception  is  in  fact  the  output  of  a  complex  system  which  uses, 
at  a  minimum,  compression  and  edge  enhancement  (32);  therefore,  added  processing  is  not  as 
unnatural  as  it  may  at  first  appear. 

To  set  the  stage  for  an  explanation  of  this  research,  a  simple  overview  of  the  human 
visual  system  is  presented  in  the  next  section. 
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2.2  Human  Visual  System 

A  diagram  of  the  human  eye  is  shown  in  figure  1.  Light  enters  through  the  cornea  and 


DIAGRAM  OF  EYE 


pupil,  gets  focused  by  the  lens,  traverses  the  vitreous  humor,  and  finally  strikes  the  retina. 
The  area  of  the  retina  containing  the  highest  concentration  of  cones  is  known  as  the  macula 
lutea,  at  the  center  of  which  is  a  pit  known  as  the  fovea  where  the  lens  of  the  eye  normally 
converges  the  incoming  light.  The  relationship  between  the  fovea  and  the  macula  lutea  is 
shown  in  figure  2.  The  retina  is  composed  of  many  layers  of  preprocessing  cells,  the  last  layer 
being  the  actual  nerve  endings.  The  multi-layer  structure  of  the  retina  is  shown  in  figure  3. 
Upon  receiving  the  signal,  the  photo-receptor  cells  send  the  signal  through  the  optic  nerve 
to  the  area  of  the  brain  responsible  for  visual  processing.  The  optical  pathway  is  shown  in 
figure  4. 

This  research  concentrates  on  enhancing  images  so  that  they  are  more  recognizable  for 
people  who  have  damage  in  the  area  of  the  macula  lutea. 
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Retina 


Figure  2.  Relationship  between  macula  lutea  and  fovea. 
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Figure  3.  Structure  of  the  retina(36) 
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Figure  4.  Optical  pathway  (39). 


2.3  Artificial  Viewing  Environment 

For  any  enhancement  techniques  to  be  practical  in  day-to-day  applications  of  face 
recognition,  a  device  must  be  employed  which  allows  the  user  to  have  a  clear  view  of  the 
algorithm  output.  The  obvious  solution  is  for  the  subject  to  carry  an  enhancement  system 
with  at  least  three  main  components:  an  input  camera,  an  enhancement  computer,  and  an 
output  screen  or  screens.  One  such  device,  the  Low  Vision  Enhancement  System  (LVES, 
pronounced  “Elvis”)  is  currently  being  employed  by  Johns-Hopkins  University.  Figure  5 
displays  this  device.  It  consists  of  three  input  cameras,  two  for  binocular  vision  and  one  for 
magnification,  and  two  miniature  LCD  output  screens  for  the  subject  to  view.  Images  are  input 
by  the  cameras  through  an  algorithm  to  be  pre-processed  in  the  hope  that  the  subject  will  be 
able  to  more  easily  discern  his  environment. 

In  his  1992  paper  concerning  limitations  of  enhancement  systems  (21),  Peli  pointed 
out  that  since  artificial  displays  have  a  finite  dynamic  range,  if  certain  frequency  components 
were  amplified  and  then  added  back  into  the  image  these  higher  intensity  values  could  exceed 
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Figure  5.  Low  Vision  Enhancement  System  (LVES)  used  as  the  input  and  output  of  an 
enhancement  system. 


the  dynamic  range  of  the  display.  This  fact,  perhaps  considered  a  limitation,  has  been  used 
to  advantage  in  some  of  Peli’s  enhancement  techniques.  If  the  image  was  simply  scaled  to 
compensate,  a  signal  amplitude  loss  occurred  at  all  frequencies. 

Also,  if  a  display  was  utilized  in  a  head-mounted  device  and  robbed  the  patient  of 
peripheral  vision,  vestibular  problems  could  occur  as  a  result  of  the  magnification  function 
associated  with  such  a  device  (21).  One  solution  to  this  could  be  to  use  a  headset  display  that 
lay  above  the  subject’s  normal  field  of  view,  one  that  could  be  referenced  at  will  and  would 
not  block  the  subjects  peripheral  cues. 

2 .4  Measurement  of  Contrast  Sensitivity 

A  subject’s  visual  contrast  threshold  at  a  specific  spatial  frequency  is  the  measure  of  the 
minimum  contrast  necessary  for  that  subject  to  perceive  a  sinusoidal  grating  at  that  frequency. 
The  contrast  sensitivity  curve  is  the  inverse  of  this  contrast  threshold  curve,  and  is  a  measure 
of  how  sensitive  or  receptive  the  subjects  visual  system  is  to  a  given  fi-equency  at  threshold 
contrast. 
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Contrast  in  digital,  gray-scale  images  is  used  to  describe  the  relationships  between  pixel 
values  in  that  image.  The  standard  measure  of  contrast  for  a  sinusoidal  grating,  an  example 
of  which  is  shown  in  figure  6,  is 


C^Lmax  Lmin^  0  <  C  <  1  (1) 

Lmax  ”1“  ^min 

where  Lmax  is  the  maximum  pixel  value  of  the  sinusoid  and  Lmin  is  the  minimum  pixel  value 
of  the  sinusoid,  usually  oscillating  about  some  mean,  DC  level. 


Figure  6.  Left:  sinusoidal  grating  with  contrast  C=l,  Right:  cross  section  of  sinusoidal 
grating. 

The  Fourier  spectrum  of  this  stimulus  was  determined  by  three  factors:  the  DC  value 
of  the  stimulus  (usually  discarded  by  the  experimenter  since  the  area  of  interest  was  in  the  AC 
frequency  spectrum),  the  frequency  of  the  sinusoid  (which  transforms  to  a  delta  function  at  the 
frequency  of  interest),  and  the  type  of  boundary  surrounding  the  sinusoid  (which  determines 
what  the  firequency  domain  delta  function  will  be  convolved  with).  If  a  sinusoidal  grating  with 
frequency  8  cycles  per  image  were  surrounded  by  a  background  with  intensity  at  the  DC  value 
of  the  grating,  and  the  sinusoid’s  boundaries  were  formed  by  straight  lines  in  the  form  of  a 
box  around  the  stimulus,  the  resulting  frequency  response  would  be  as  shown  in  figure  7.  This 
indicates  that  the  measurement  of  the  response  of  the  subject  to  this  stimulus  sheds  light  on  the 
processing  of  the  non-zero  areas  in  the  graph  on  the  right  in  figure  7.  To  concentrate  the  area 
of  interest  more  to  the  frequencies  around  8  cycles  per  image,  a  sinusoid  modulated  by  a  two 
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Cycles  per  image 


Cycles  per  Image 


0  0 

Figure?.  Left:  boxed  sinusoidal  grating  with  contrast  C=l.  Right:  closeup  of  frequency 
spectrum  to  include  the  area  of  interest,  here  8  cycles  per  image. 

dimensional  Gaussian  rather  than  a  two  dimensional  rect  function  is  sometimes  used  (25,  3). 
The  result  of  this  operation  is  shown  in  figure  8. 

2.4.1  The  staircase-method  in  psychophysics.  When  looking  for  a  psychometric 
threshold,  such  as  the  minimum  amplitude  of  a  sinusoidal  grating  at  a  specific  spatial  frequency 
so  that  the  grating  is  seen  by  a  subject,  Comsweet  argued  that  the  staircase-method  provides  a 
powerful  and  efficient  estimation  tool  (6).  For  the  above  example,  a  subject  would  be  presented 
with  a  grating  and  respond  with  a  positive  or  negative  as  to  whether  the  stimulus  was  seen. 
If  the  answer  was  positive,  the  contrast  of  the  sinusoid  would  be  lowered;  if  the  answer  was 
negative,  the  contrast  would  be  increased.  In  the  ideal  case,  a  subject’s  responses  would 
simply  oscillate  about  the  threshold  amplitude  as  shown  in  figure  9,  and  the  threshold  level 
would  be  estimated  to  lie  between  the  points  of  the  first  reversal.  Comsweet  points  out  that 
in  actual  experiments  a  subject’s  responses  were  not  independent  of  one  another  and  results 
obtained  when  trying  to  determine  a  visual  threshold  using  the  staircase  method  would  more 
closely  match  those  shown  in  figure  10.  He  said  the  actual  threshold  was  usually  estimated  as 
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Figure  8.  Left:  sinusoidal  grating  with  contrast  C=1  modulated  by  a  two  dimensional  Gaus¬ 
sian  function.  Right:  closeup  of  frequency  spectrum  to  include  the  area  of  interest, 
once  again  8  cycles  per  image. 

the  intensity  level  above  which  50%  of  the  responses  were  “yes”  and  below  which  50%  of  the 
responses  were  “no”. 

The  experimenter  had  control  over  four  parameters  when  using  the  staircase-method: 
the  initial  amplitude  or  contrast,  the  step-sizes  between  subsequent  trials,  how  many  trials 
until  the  test  was  terminated,  and  modifications  of  step-size  if  warranted.  The  first  few 
trials  of  a  test  such  as  this  were  necessary  to  establish  the  approximate  range  in  which  the 
threshold  lay;  therefore,  the  initial  amplitude  or  contrast  was  chosen  as  closely  as  possible  to 
the  estimated  threshold  so  that  as  few  trials  were  made  as  necessary.  The  step-size  between 
trials  was  ideally  chosen  to  approximate  sensory  perception  intervals  that  occur  in  nature. 
How  many  trials  used  depended  on  the  experimenter’s  preference  since  a  dichotomy  existed 
between  expediency  (low  number  of  trials)  and  accuracy  (high  number  of  trials).  As  the  test 
progressed,  the  experimenter  could  alter  step  sizes  in  a  heuristic  fashion  (see  the  next  section 
concerning  parameter  estimation  by  sequential  testing)  to  home  in  on  the  actual  threshold 
being  estimated. 


Figure  9. 


Figure  10. 
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•  "NO" 

Trial  Number  x  'yes" 


Ideal  results  of  using  the  staircase-method  to  determine  visual  contrast  thresh¬ 
old  (6). 


Typical  results  of  using  the  staircase-method  to  determine  visual  contrast  thresh¬ 
old  (6). 
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One  of  the  problems  with  this  method  as  pointed  out  by  Comsweet  was  that  if  the 
answers  at  each  trial  consisted  of  a  simple  yes/no  answer,  a  subject  could  essentially  “close 
their  eyes”  as  the  testing  proceeded  and  “game”  the  test  (whether  or  not  they  realized  they 
were  doing  it)  so  that  the  results  would  be  unreliable.  He  conjectured  that  the  double-staircase 
method  offered  a  way  to  circumvent  this  problem.  The  first  trial  presented  a  testing  level  on 
one  side  of  the  estimated  threshold  while  the  second  trial  presented  a  testing  level  on  the  other 
side.  See  figure  11.  Every  other  trial,  then,  belonged  to  a  distinct  set  of  two  trials  occurring 
simultaneously.  To  add  another  dimension,  the  test  could  randomly  decide  to  which  set  a  trial 
belonged  as  the  test  commenced  as  in  figure  12.  Comsweet  argued  that  the  subject  must  then 
be  responding  to  the  actual  stimuli  and  not  guessing. 
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Figure  11.  Possible  results  of  using  the  double  staircase  method  to  determine  threshold  (6). 


2.4.2  Parameter  Estimation  by  Sequential  Testing  (PEST).  Taylor  expanded  upon 
the  staircase-method  by  providing  a  set  of  mles  based  on  the  Wald  sequential  likelihood-ratio 
test  (34)  as  to  when  the  testing  level  should  be  changed  and  the  magnitude  of  that  change  (33). 
The  PEST  algorithm  was  to  be  an  adaptive  method  of  testing  for  a  psychometric  function  in  that 
the  step  size  and  number  of  trials  would  not  depend  exclusively  on  predetermined  parameters, 
but  would  instead  depend  on  previous  trials  in  the  same  experiment.  Taylor  pointed  out  that 
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E[N{C)]  =PtXT±W.  (2) 

N(C)  was  the  total  number  of  correct  responses  after  each  change  in  testing  level,  E[...]  was 
the  expected  value  operator,  and  W  was  a  constant  determined  heuristically  based  on  how 
powerful  the  test  was  to  be.  The  deviation  limit,  W,  determined  the  power  versus  speed 
dichotomy  of  the  algorithm.  Large  W  resulted  in  many  trials  and  estimates  very  close  to  the 
actual  threshold  level,  while  small  W  caused  rapid  testing  at  the  expense  of  accuracy.  Pt 
was  the  target  probability;  this  target  probability  was  the  ratio  of  the  estimated  number  of 
correct  responses  to  the  number  of  trials  at  that  particular  testing  level.  For  example,  the 
experimenter  could  say  he  felt  getting  5  out  of  6  responses  correct  would  indicate  that  the 
subject  could,  in  fact,  discriminate  a  horizontal  grating  from  a  vertical  grating  in  a  measure 
of  contrast  sensitivity.  He  would  then  set  Pt  =  5/6  or  .8333.  The  value  chosen  for  Pt  did  not 
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effect  the  validity  of  the  final  estimate  of  the  threshold  level,  but  it  did  effect  the  number  of 
trials  necessary  to  arrive  at  that  estimate. 

In  the  context  of  testing  for  a  contrast  threshold  curve  for  a  given  frequency,  grating 
contrast  of  the  next  trial  was  based  on  preceding  trials.  Taylor’s  rules  were  outlined  as  follows: 

•  When  the  above  boundary  rule  indicated  a  reversal  of  direction;  for  example,  the 
contrast  stopped  going  up  and  on  the  next  step  went  down;  the  contrast  step  size  was  halved. 

•  When  the  boundary  rule  indicated  the  second  step  to  be  in  the  same  direction  as  the 
first,  the  second  step  size  was  the  same  as  the  first. 

•  The  size  of  the  fourth,  fifth,  sixth,  and  so  on,  step  in  the  same  direction  were  double 
the  preceding  step  size. 

•  The  step  size  of  the  third  successive  step  in  the  same  direction  depended  on  the  last 
reversal  preceding  the  current  step  direction  series;  if  the  step  size  preceding  that  reversal  was 
double  its  preceding  step,  this  step  size  was  the  same  as  the  second  step  size;  if  the  step  size 
preceding  that  reversal  was  not  double  its  preceding  step,  this  step  size  was  double  the  second 
step  size. 

The  PEST  terminated  when  a  step  was  called  for  that  equalled  the  experimenter ’s  chosen 
minimum  step  size.  The  estimated  value  of  the  contrast  threshold  level  being  sought  was  then 
the  next  level  called  for  after  the  last  test.  No  testing  occurred  at  that  level.  Methods  have 
since  been  outlined  by  Pentland  (28,  16)  which  minimize  the  number  of  trials  necessary  to 
obtain  a  given  threshold  measurement  using  PEST. 

Peli  generated  contrast  threshold  curves  for  several  subjects  (19).  These  are  shown  in 
figure  13  and  were  generated  using  the  PEST  algorithm  (33, 28)  and  sinusoidal  Gabor  patches 
(as  shown  in  figure  8)  at  frequencies  of  2, 4,  8, 16,  and  32  cycles  per  image  (cpi). 

2.5  Measurement  of  Suprathreshold  Characteristics 

Cannon  demonstrated  that  visual  perception  contrast  curves  generated  by  measuring 
contrast  above  threshold  (suprathreshold)  characteristics  of  an  observer  were  not  simply 
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Figure  13.  Various  measured  threshold  curves  of  persons,  the  curve  used  for  simulation 
purposes  was  the  one  corresponding  to  20/136  vision  (24) 


contrast  threshold  curves  shifted  upward  (2).  He  showed  that  as  the  overall  contrast  was 
increased  for  all  spatial-frequency  sinusoidal  gratings,  the  amount  of  contrast  necessary  to 
make  two  spatial  frequencies  look  as  though  they  had  the  same  contrast  approached  the  same 
value.  Figure  14  displays  what  these  curves  look  like  in  normal  individuals.  For  example; 
using  one  of  the  curves  from  figure  14,  a  sinusoidal  grating  varying  at  4  cycles  per  deg  (cpd) 
was  displayed  at  a  contrast  of  0.010  and  another  grating  varying  at  16  cycles  per  deg  was 
displayed  simultaneously  at  some  initial,  different,  contrast.  The  16  cpd  grating  was  adjusted 
by  the  subject  until  that  subject  perceived  the  two  contrasts  to  be  equal.  On  this  particular 
curve,  the  actual  contrast  of  the  16  cpd  grating  was  around  0.020.  The  subject  perceived  the 
two  contrasts  to  be  equal  even  though  the  actual  contrast  of  the  16  cpd  grating  was  twice  as 
high  as  the  4  cpd  grating.  At  threshold,  the  subject  perceived  the  4  cpd  grating  at  an  actual 
threshold  contrast  and  the  16  cpd  grating  at  a  different  actual  threshold  contrast.  As  the  4  cpd 
grating  contrast  was  increased,  the  perceptible  difference  between  the  16  cpd  actual  contrast 
and  the  4  cpd  actual  contrast  decreased.  This  eventually  led  to  a  minimum  contrast  at  which 
the  actual  contrasts  of  different  spatial  frequencies  with  the  same  contrast  were  perceived  to 
be  the  same. 
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Figure  14.  Example  of  suprathreshold  curves  in  individuals  with  normal  vision  (2) 
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The  method  employed  by  Cannon  to  measure  contrast  suprathreshold  characteristics  of 
a  subject  was  more  straightforward  than  was  the  measurement  of  contrast  threshold  charac¬ 
teristics.  The  subject  was  presented  a  grating  at  a  reference  spatial  frequency.  Other  spatial 
frequencies  were  presented  and  contrast  matched  by  the  subject  to  this  reference  frequency. 
The  contrast  of  the  reference  frequency  was  increased  for  each  set  of  tests  until  a  flat  curve 
was  obtained. 

2.6  Simulation  of  What  a  Person  With  a  Given  Threshold  Curve  Sees  When  Looking  at  an 
Image 

To  simulate  how  the  human  visual  system  processes  information  contained  in  an  image, 
Peli  proposed  a  nonlinear  model  using  the  concept  of  local  band-limited  contrast  (19).  This 
model  generated  a  simulation  of  what  a  person  with  a  given  threshold  curve  saw  when  looking 
at  an  image.  Figure  15  displays  the  128x128  image  presented  to  the  observer  in  this  example 


Figure  15.  Original  Image  (size  128  x  128) 

simulation.  The  dimensions  of  the  image  in  this  figure  are  1 .257  inches  square.  If  the  reader 
holds  the  figure  1.5  feet  from  the  eyes,  the  image  subtends  4  degrees  of  arc;  one  cpd  of  the 
image  is  four  cpi. 

The  first  step  in  the  simulation  was  to  break  the  image  into  its  constituent  band-pass 
filtered  components,  much  as  Ginsburg  conjectured  the  brain  does  (8).  There  were  four  images 
associated  with  each  of  the  five  spatial  frequency  contrast  threshold  responses  obtained  with 
the  contrast  sensitivity  testing:  four  images  associated  with  2  cpi,  four  with  4  cpi,  four  with 
8  cpi,  four  with  16  cpi,  and  four  with  32  cpi.  The  first  image  generated  was  obtained  by 
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filtering  the  original  image’s  frequency  spectrum  with  a  concentric  log-cosine  filter  of  the 
form: 


Gi{r)  =  1/2[1  -I-  cos{'Klog2r  —  i  =  0, 1, 2, ...  (3) 


Plots  of  the  band-pass  filters  centered  at  2, 4, 8, 16,  and  32  cpi  are  shown  on  the  left  in  figure  16, 
and  the  related  low-pass  filters  designed  to  capture  the  energy  under  their  corresponding  band¬ 
pass  filter  are  shown  on  the  right  of  that  same  figure.  The  text  in  the  plot  indicates  where  the 
center  frequency  of  the  band-pass  filter  lay  and  where  the  related  low-pass  filter  became  zero. 
Each  band-pass  filter  was  centered  at  2*  cpi  and  was  one  octave  wide.  Using  16  cpi  as  an 
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Figure  16.  Cross-section  log-cosine  filters  necessary  for  Peli’s  simulation. 

example,  the  filter  G4  was  centered  at  16  cpi  and  spanned  from  8  to  32  cpi,  with  the  amplitude 
of  the  filter  at  half  maximum  at  frequencies  of  12  and  24  cpi  (thus,  the  octave  bandwidth). 
For  the  simulation,  the  image  generated  by  band  pass  filtering  (center  frequency  16  cpi)  the 
original  image  is  shown  on  the  left  in  figure  17  and  was  point-by-point  divided  by  the  image 
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shown  on  the  right  in  figure  17  which  was  generated  by  low-pass  filtering  the  original  image 
by  the  corresponding  low-pass  filter.  The  resulting  image  from  this  division  was  defined  by 


Figure  17.  Left:  image  band-pass  filtered  with  center  of  filter  at  16  cycles  per  image  (DC 
value  of  128  added  for  image  perception).  Right:  corresponding  low-pass-filtered 
image. 

Peli  as  a  “local  band-limited  contrast  image,”  and  each  pixel  value  in  each  image  was  in  fact  a 
contrast  in  that  band  (in  this  case,  16  cpi)  at  that  location.  A  local  band-limited  contrast  image 
was  generated  for  2,  4,  8,  16,  and  32  cycles  per  image.  The  pixel  values,  or  contrasts,  in  each 
contrast  image  were  compared  to  the  threshold  value  measured  in  the  contrast  sensitivity  test 
corresponding  to  the  appropriate  frequency.  For  example,  if  the  contrast  threshold  of  a  subject 
at  16  cycles  per  image  was  .55,  then  each  pixel  value  in  the  contrast  image  corresponding  to  16 
cycles  per  image,  as  in  the  above  example,  were  compared  to  .55.  If  the  pixel  value  was  less 
than  .55,  the  pixel  value  at  that  location  in  the  image  generated  by  the  16  cpi  band-pass  filter 
(left  image  in  figure  17)  was  set  to  zero  since  the  contrast  at  that  location  in  that  frequency 
band  was  below  the  observable  threshold  of  the  subject.  If  the  pixel  value  was  greater  than  or 
equal  to  .55,  the  value  at  that  location  in  the  bpf  image  was  left  alone.  The  images  obtained 
by  this  process  were  called  threshold  images,  and  were  added  together  along  with  the  lowest 
frequency  component  to  complete  the  simulation.  Peli  reasoned  that  this  resulting  image 
was  representative  of  what  a  person  with  that  particular  threshold  response  would  see  when 
looking  at  figure  15.  The  threshold  values  used  by  Peli  for  his  initial  simulations  were  obtained 
from  a  subject  with  20/136  vision  shown  on  the  graph  in  figure  13.  Here,  that  same  curve 
was  used;  the  image  produced  is  shown  in  figure  18  and  represented  what  that  low-vision 
observer  perceived  when  looking  at  the  original  image.  The  Matlab  code  used  to  accomplish 
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this  simulation  is  shown  in  Appendix  A.  To  assure  that  the  simulation  was  accurate  for  a 
person  with  normal  vision,  the  curve  representing  normal  vision  in  figure  13  was  used  in 
the  simulation.  The  resulting  picture  was  as  shown  in  figure  19.  Notice  that  while  the  high 
frequencies  have  been  attenuated,  the  face  is  still  recognizable  (assuming,  of  course,  that  the 
reader  has  “normal”  vision). 


Figure  1 8 .  The  image  on  the  right  is  a  simulation  of  what  a  person  with  the  curve  correspond¬ 
ing  to  20/136  visual  acuity  shown  in  figure  13  on  page  18  sees  when  looking  at 
image  on  left. 


Figure  19.  The  image  on  the  right  is  a  simulation  of  what  a  person  with  the  curve  corre¬ 
sponding  normal  visual  acuity  shown  in  figure  13  on  page  18  sees  when  looking 
at  image  on  left. 


2.7  Digital  Enhancement 

Peli  and  Peli  first  introduced  the  concept  of  using  digital  image  processing  to  enhance 
images  presented  to  patients  suffering  from  visual  degradation  in  1984  (26).  The  initial  model 
of  a  person’s  visual  loss  they  used  was  the  contrast  sensitivity  function  (CSF)  which  was  the 
inverse  of  the  contrast  threshold  curve  which  represented  the  minimum  contrast,  or  amplitude, 
necessary  for  a  person  to  resolve  whether  sinusoidal  intensity  gratings  at  various  frequencies 
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were  horizontally  or  vertically  oriented.  They  measured  the  CSF  in  various  patients  and  used 
them  to  calculate  individual  visual  degradation  transfer  functions  (VDTF),  which  were  the  ratio 
of  the  CSF  of  the  patients  with  cataractous  vision  to  the  CSF  of  persons  with  normal  vision. 
Peli  decided  that  using  the  VDTF  as  an  inverse  filter  produced  too  much  “salt  and  pepper” 
noise  due  to  the  amplification  of  the  higher  frequencies  attenuated  by  a  degraded  visual  system. 
He  instead  conjectured  that  adaptive  image  enhancement  (27)  (see  below)  may  be  better  for 
enhancement.  Although  no  testing  was  done  in  1984  to  see  the  effects  of  these  enhancements 
on  patients,  the  effects  of  inverse  filter  enhancement  as  viewed  through  simulated  cataractous 
conditions  were  studied.  In  their  paper  on  computerized  image  enhancement  (22),  Peli 
and  others  addressed:  contrast  manipulation,  spatial  filtering,  and  pseudocolor  recoding.  In 
her  paper  on  attempts  to  improve  reading  rates  in  the  visually  impaired,  Lawton  discussed 
using  the  contrast  sensitivity  function  as  a  compensation  filter  (15).  Peli  once  again  tackled 
enhancement  using  a  pre-emphasis  filter  concept  (2 1 )  in  his  1 992  paper  on  limitations  of  image 
enhancement. 

2.7.1  Contrast  Manipulation.  Contrast  manipulation  involved  using  the  dynamic 
range  of  a  display  device  to  redistribute  the  contrast  information  in  such  a  way  so  that  it  was 
more  useful  to  the  observer.  Two  types  of  contrast  manipulation  were  “histogram  equalization 
and  hyperbolization”,  and  “thresholding”.  Histogram  equalization  forced  the  image’s  gray 
scale  values  to  yield  a  flat  histogram.  The  image  on  the  right  in  figure  20  is  the  histogram 
equalized  version  of  the  low  contrast  image  on  the  left  in  the  same  figure.  Peli  pointed  out, 
though,  that  histogram  equalization  could  sometimes  cause  a  loss  of  information  since  bright 


Figure  20.  Result  of  histogram  equalization  on  a  low  contrast  image. 
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regions  in  the  original  image  may  become  very  bright  after  equalization.  Hyperbolization 
forced  the  values  into  a  hyperbolic  shaped  histogram  and  had  less  of  the  brightening  effect  of 
histogram  equalization.  Thresholding,  used  in  many  applications  for  segmenting  information 
in  an  image,  yielded  the  outline  of  the  person  in  the  original  image  and  preserved  much  of 
the  information  as  shown  in  the  left  image  in  figure  21.  The  middle  image  of  figure  21  is 
the  result  of  passing  the  image  through  the  previous  simulation  of  what  a  person  with  20/136 
vision  sees  when  looking  at  the  image  on  the  left.  Adaptive  thresholding  (24)  involved 


Figure  21.  Left:  Result  of  thresholding  in  original  image.  Middle:  Simulation  of  what  a 
person  with  the  curve  corresponding  to  20/136  visual  acuity  shown  in  figure  13 
on  page  18  sees  when  looking  at  image  on  left.  Right:  Simulation  of  what  that 
same  person  sees  when  looking  at  the  original,  un-enhanced  image. 

dividing  an  image  into  non-overlapping  sections.  Each  section  was  evaluated  to  determine  its 
appropriate  threshold  and  every  pixel  in  the  image  was  then  thresholded  based  on  threshold 
values  interpolated  between  the  sections.  Peli  stated  one  drawback  to  this  method  was  that 
some  patients  complained  that  thresholding  as  an  enhancement  technique  caused  the  resulting 
image  to  look  “cartoonish.” 

2.7.2  Spatial  Filtering.  Image  enhancement  by  spatial  filtering  would  optimally 
amplify  those  spatial  frequencies  most  attenuated  by  the  visual  system  of  the  afflicted  patient. 
At  the  time  of  this  thesis,  Peli  was  investigating  the  use  of  the  VDTF  in  an  inverse  filter 
enhancement  technique  (21).  He  mentioned  a  number  of  enhancements  which  would  boost 
the  high  frequency  information  in  an  image.  High-pass  image  filtering  techniques  included: 
unsharp  masking,  extremum  sharpening,  and  adaptive  filtering.  Unsharp  masking  subtracted 
a  low-pass  version  of  the  image  from  the  original.  Extremum  sharpening  passed  a  spatial 
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mask  over  the  original  image;  the  points  inside  this  mask  were  assigned  the  maximum  or 
minimum  pixel  value  in  the  window,  depending  on  which  value  the  current  pixel  was  closest 
to.  Adaptive  irnage  enhancement  involved  low-pass  filtering  the  original  image,  subtracting 
this  Ipf  image  from  the  original  image,  and  adding  this  new  image  (which  contained  high 
frequency  information)  to  the  original.  The  effects  of  this  technique  are  shown  in  figure  22; 
on  the  left  of  this  figure  is  the  enhanced  image,  while  in  the  middle  is  the  enhanced  image  as 


Figure  22.  Left:  Result  of  Adaptive  filtering  original  image.  Middle:  Simulation  of  what  a 
person  with  the  curve  corresponding  to  20/136  visual  acuity  shown  in  figure  13 
on  page  18  sees  when  looking  at  image  on  left.  Right:  Simulation  of  what  that 
same  person  sees  when  looking  at  the  original,  unenhanced  image. 

seen  by  the  low  vision  observer  with  visual  acuity  of  20/136  shown  in  figure  13.  Comparing 
the  image  in  the  middle  to  the  unenhanced  perceived  image  on  the  right,  notice  that  many 
of  the  facial  characteristics  have  been  recovered  and  that  the  enhanced  image  may  be  more 
recognizable  to  the  low  vision  patient.  Adaptive  filtering  was  the  same  as  adaptive  image 
enhancement. 

Peli  tested  adaptive  image  enhancement  and  adaptive  thresholding  in  1991  on  patients 
suffering  from  low  vision  (24).  He  used  parametric  estimation  by  sequential  testing  (PEST) 
(28,  16)  to  minimize  the  number  of  measurements  necessary  to  determine  each  patient’s  CSF. 
These  enhancement  techniques  were  not  tuned  for  a  specific  patient’s  CSF,  but  were  tuned  to 
enhance  the  frequencies  above  4  cycles  per  image  which  were  most  useful  for  face  recognition. 
The  results  of  the  testing  indicated  that  39  out  of  46  patients  tested  had  better  facial  recognition 
ability  when  viewing  the  enhanced  faces  than  when  viewing  the  un-enhanced  faces. 
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2.7.3  Pseudocolor  Recoding.  Images  could  be  enhanced  using  pseudocolor  re¬ 
coding,  which  would  take  advantage  of  the  extra  information  contained  in  the  “hues”  of  the 
image.  Subjects,  though,  would  need  to  be  retrained  to  the  new  color  representations.  Also, 
damage  of  the  macula  lutea  implies  damage  near  the  fovea,  the  area  containing  the  highest 
concentration  of  the  cones  which  detect  color.  Therefore,  this  technique  would  probably  prove 
itself  less  than  optimal  for  improving  recognition  ability  in  those  suffering  from  low-vision 
due  to  macular  degeneration. 

2.7.4  Compensation  Filters.  In  her  1988  paper  (15),  Lawton  described  her  attempt 
at  improving  word  recognition  in  three  low-vision  patients  using  compensation  filters  based 
directly  on  each  individual  patient’s  CSF  normalized  to  a  normal  viewer’s  CSF.  She  used 
this  normalized  CSF  (NCSF)  essentially  as  an  inverse  filter,  directly  correlating  the  amount 
of  gain  at  a  particular  frequency  to  the  inverse  of  the  NCSF  at  that  spatial  frequency.  The 
compensation  filter  used  was  of  the  form. 


NCSF 

ncsWTTS’ 

where  Gain  determined  how  much  of  a  tradeoff  there  was  to  be  between  amplification  of  high 
frequency  noise  and  a  noise-free  image.  In  her  experiment,  only  letters  were  enhanced,  not 
faces.  Testing  of  this  technique  revealed  that  enhancing  the  words  with  the  compensation  filters 
allowed  recognition  that  would  have  required  a  27-70%  magnification  of  the  un-enhanced 
word. 

The  left  image  in  the  figure  23  was  obtained  by  using  a  compensation  filter  generated 
from  the  two  curves  in  figure  1 3  on  page  1 8,  one  representing  normal  vision  and  one  represent¬ 
ing  the  vision  of  the  person  with  20/136  visual  acuity.  The  NCSF  was  the  normal  threshold 
curve  point-by-point  divided  by  the  low-vision  curve.  The  gain  in  equation  4  was  set  at  1000, 
which  caused  the  compensation  filter  to  behave  as  an  inverse  filter.  To  avoid  divide-by-zero 


27 


Figure  23.  Left:  Result  of  applying  compensation  filter  on  original  image.  Middle:  Simula¬ 
tion  of  what  a  person  with  the  curve  corresponding  to  20/136  visual  acuity  shown 
in  figure  13  on  page  18  sees  when  looking  at  image  on  left.  Right:  Simulation  of 
what  that  same  person  sees  when  looking  at  the  original,  unenhanced  image. 


problems,  frequencies  above  32  cpi  were  set  to  zero  during  the  filter  generating  process.  A 
cross  section  of  the  compensation  filter  used  for  enhancement  is  shown  in  figure  24.  The 


Figure  24.  Cross  section  of  compensation  filter  used  for  enhancement. 

resulting  enhanced  image  was  histogram  equalized  so  that  the  pixel  values  would  lie  between 
0  and  500.  This  resulted  in  saturation  of  some  of  the  pixels  since  the  picture  was  displayed 
with  a  dynamic  range  between  0  and  255. 

In  1994,  Fine  and  others  attempted  to  reproduce  the  results  obtained  by  Lawton  in  her 
experiment  concerning  improved  reading  rates  in  the  visually  impaired  (7).  Although  they 
could  not  reproduce  her  results  precisely,  they  were  able  to  find  a  significant  increase  in  the 


reading  rates  of  individuals  viewing  enhanced  text .  Fine  used  letters  similar  in  appearance  to 
Lawton’s.  Conjecturing  that  Lawton’s  enhancement  technique  may  have  enlarged  the  letters, 
they  set  up  a  test  to  find  out  if  larger  letters  would  lead  to  the  same  increase  in  reading 
rate  as  Lawton’s  technique.  It  was  seen  that  simply  enlarging  the  letters  did  not  lead  to  the 
same  increase  in  reading  rate,  thereby  arguing  the  possible  usefulness  of  using  the  normalized 
contrast  sensitivity  function  as  an  inverse  filter. 

2.7.5  Modified  Pre-Emphasis  Model.  Peli’s  most  recent  image  enhancement 
endeavor  involved  using  a  band-emphasis  filter  to  boost  frequencies  most  useful  for  face 
recognition  in  a  given  subject  (21).  The  image  was  first  low-pass  filtered  to  remove  all 
frequencies  unusable  by  the  subject.  An  emphasis  filter  was  built  to  boost  the  highest  frequency 
seen  by  the  subject,  while  leaving  untouched  the  frequencies  below  that  highest  frequency. 
The  amount  to  boost  the  high  frequency  was  determined  by  the  ratio  of  the  subject’s  contrast 
threshold  to  a  normal  contrast  threshold  at  that  frequency.  Peli  left  the  lower  frequencies 
untouched  for  aesthetic  reasons  since  those  frequencies  kept  the  enhanced  image  from  looking 
like  a  line-drawing.  The  filter  used  is  shown  in  figure  25.  It  was  also  pointed  out  that  some 


Figure  25.  Pre-emphasis  filter  used  for  enhancement. 

saturation  of  the  pixel  values  in  the  resulting  enhanced  image  allowed  more  information  to 
be  displayed  than  if  straight  histogram  equalization  to  the  full  dynamic  range  of  the  display 
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was  used.  As  an  example,  if  the  frequencies  of  the  original  image  are  multiplied  by  the 
pre-emphasis  filter  shown  in  figure  25  and  then  histogram  equalized  to  stretch  the  pixel  values 
between  0  and  800,  the  left  image  in  figure  26  is  the  result.  That  image  is  the  enhanced  image 


Figure  26.  Left:  Result  of  applying  pre-emphasis  filter  on  original  image.  Middle:  Simula¬ 
tion  of  what  a  person  with  the  curve  corresponding  to  20/136  visual  acuity  shown 
in  figure  13  on  page  18  sees  when  looking  at  image  on  left.  Right:  Simulation  of 
what  that  same  person  sees  when  looking  at  the  original,  unenhanced  image. 

presented  to  the  subject,  while  the  image  in  the  middle  is  the  simulated  image  perceived  by  a 
person  with  20/136  vision  shown  in  the  threshold  curve  of  figure  13.  Once  again,  notice  that 
the  middle  image  may  be  more  recognizable  than  the  un-enhanced  perceived  image  on  the 
right. 

2.8  Wavelets 

The  previously  mentioned  enhancement  techniques  represent  a  large  and  useful  arsenal 
in  the  fight  to  restore  recognition  abilities  in  those  with  low-vision.  The  most  promising  of  the 
above  rely  on  direct  manipulation  of  an  image’s  spatial  frequency  spectrum.  Some  theories 
involving  the  way  mammalian  brains  process  auditory  (35)  and  visual  (12,  13)  information 
suggest  that  human  perception  can  be  modelled  with  wavelet  analysis.  The  significance  of 
discrete  wavelet  analysis  compared  to  discrete  Fourier  analysis  lies  in  the  decomposition  and 
reconstruction  process  inherent  to  wavelet  analysis.  Here,  it  is  seen  that  discrete  wavelet 
analysis  can  be  used  as  a  tool  for  image  enhancement  and  for  simulating  human  perception  of 
an  image  using  a  particular  contrast  sensitivity  function. 
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In  his  paper  on  contrast  in  complex  images,  Peli  mentioned  image  analysis  and  syn¬ 
thesis  schemes  using  his  log  cosine  implementation  to  determine  band-limited  contrast  (24). 
One  analysis  scheme  that  has  not  been  used,  one  that  may  more  closely  represent  visual 
processing,  is  wavelet  analysis.  Discrete  wavelet  coefficient  sets  can  be  thought  of  as  images 
that  have  been  filtered  in  the  frequency  domain  with  filters  having  cutoffs  at  powers  of  two 
and  bandwidths  of  one  octave.  This  section  explains  discrete  wavelet  implementation  from  a 
non-mathematically  rigorous,  very  practical  standpoint.  For  a  more  mathematically  intense 
analysis  of  two  dimensional  wavelet  theory,  an  excellent  overview  is  presented  in  the  paper 
by  Mallat  (18)  and  in  the  Masters  Thesis  by  Laing  (14)  (A  note  concerning  the  Laing  Thesis; 
in  figures  36  and  38,  both  explaining  the  mechanics  of  wavelet  reconstruction,  the  upsampling 
occurs  after  convolution.  This  is  incorrect,  although  it  is  explained  correctly  in  his  text;  the 
convolution  must  occur  after  the  upsampling,  thus  maintaining  the  reconstruction  as  the  mirror 
operation  of  the  decomposition). 

To  decompose  an  image  using  the  standard  discrete  wavelet  transform  (DWT),  two 
one-dimensional  vectors  are  necessary:  a  vector  representing  the  scaling  function  and  a 
vector  representing  the  wavelet  function.  The  vector  representing  the  scaling  function  is 
usually  denoted  as  h  and  has  a  frequency  response  displaying  low-pass  characteristics.  The 
vector  representing  the  wavelet  function  is  usually  denoted  as  g  and  has  a  frequency  response 
displaying  high-pass  characteristics,  h  and  g  form  an  orthogonal  basis  set. 

With  simple  orthogonal  wavelets,  a  single  vector  is  given  for  the  scaling  function.  The 
wavelet  coefficients  vector  is  generated  from  the  scaling  coefficients  using  the  relationship. 


g[n)  =  —1^  ”  X  h{l  —  n).  (5) 

Now  that  these  two  vectors  are  available,  the  actual  image  decomposition  is  accomplished  by 
using  the  mirrors  of  these  two  vectors;  that  is,  h{x)  =  h{—x)  and  g{x)  =  g{—x)  are  the 
vectors  that  are  actually  used  in  the  decomposition  operation  shown  in  figure  28. 
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The  wavelet  set  used  in  this  research  was  the  biorthogonal  wavelet  set  used  by  the 
FBI  for  image  compression  of  fingerprint  data  (10,  37,  17).  This  set  was  chosen  because  the 
decomposition  functions  were  symmetric  about  zero,  which  meant  that  manipulations  of  the 
detail  coefficients  would  not  cause  distortion  of  the  reconstructed  image  (as  may  be  expe¬ 
rienced  using  other  compactly  supported,  non-symmetric  wavelet  sets).  Using  biorthogonal 
wavelets  required  a  slightly  different  approach  than  that  for  orthogonal  wavelets  sets  outlined 
in  Laing’s  thesis.  With  a  biorthogonal  wavelet  set,  two  vectors  (used  for  decomposition)  are 
initially  given.  One  vector  given  represents  the  scaling  function  to  use  for  decomposing  an 
image  (h)  and  the  other  vector  given  represents  the  wavelet  function  used  for  decomposition 
(g).  The  biorthogonal  vectors  shown  in  tables  1  and  2  were  the  ones  used  by  the  FBI  for 
fingerprint  encoding  and  were  obtained  from  the  Thesis  by  MacDonald  (17).  The  waveforms 
corresponding  to  these  vectors,  used  for  decomposition,  and  their  relationship  are  shown  in 
figure  27.  Notice  that  figure  27  seems  to  indicate  that  non-zero  values  for  the  scaling  and 

Table  1.  Coefficients  for  the  scaling  function  used  in  the  forward  wavelet  decomposition  (h) 


Position 

Approximate  Value 

ho 

-1-0.85269867900940 

h±i 

-1-0.37740285561265 

h±2 

-0.11062440441842 

h±3 

-0.023849465019380 

h±i 

-t-0.037828455506995 

Table  2.  Coefficients  for  the  wavelet  function  used  in  the  forward  wavelet  decomposition  (g) 


Position 

Approximate  Value 

9-1 

-1-0.78848561640566 

9-2,0 

-0.41809227322221 

9-3,1 

-0.040689417609558 

9-4,2 

-1-0.064538882628938 

wavelet  functions  occur  only  at  positive  indices.  This  is  an  artifact  of  using  Matlab  for  the 
analysis.  The  scaling  function  (h)  is  symmetric  about  zero,  and  the  wavelet  function  (g)  is 
symmetric  about  minus  one.  Figure  28  shows  the  decomposition  process  for  an  image. 

Decomposing  an  N  x  N  image  involved  implementing  the  following  steps: 
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Figure  27.  Decomposition  scale  function  (h)  is  on  the  left,  decomposition  wavelet  function 
(g)  is  in  the  middle,  while  their  relationship  is  shown  on  the  right. 

1.  Generate  a  column  coefficient  set  by  convolving  the  columns  of  the  image  with 
generate  another  column  coefficient  set  by  convolving  the  columns  of  the  image  with 
h.  The  convolution  of  the  first  column  in  the  image  with  the  scale  function  is  defined  as 

_  AT 

Columui  *h  =  E  C  oluTflTluih^—jYi,  (6) 

m=l 

This  convolution  is  circular;  that  is  to  say  that  when  the  vector  encounters  the  image’s 
edge,  it  wraps  around  to  the  other  side  of  the  image. 

2.  Downsample  the  columns  of  each  of  these  coefficient  sets  by  keeping  every  other  row. 
The  size  of  each  image  matrix  is  now  A^/2  x  N. 

3.  Generate  two  coefficient  sets  by  convolving  the  rows  of  each  of  the  two  downsampled 
column  coefficient  sets  with  g;  generate  two  more  coefficient  sets  by  convolving  the 
rows  of  each  of  the  two  downsampled  column  coefficient  sets  with  h. 

4.  Downsample  the  rows  of  each  of  these  coefficient  sets  by  keeping  every  other  column. 
The  sizes  of  each  of  the  four  new  matrices  is  now  N/ 2  x  N/ 2. 

The  detail  coefficients  obtained  during  a  wavelet  decomposition  can  be  thought  of  in 
the  frequency  domain  as  multiplying  the  high  frequency  components  of  the  original  image  by 
a  one-octave-wide  band-pass  filter  (5).  The  cutoffs  of  this  bpf  are  always  a  power  of  two  when 
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:  Convolve  with  response  function  xfnj 


Figure  28.  Decomposition  flow  chart  to  obtain  first  level  of  decomposition 
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the  original  matrix  has  dimensions  that  are  powers  of  two.  For  example,  the  image  shown 
in  figure  15  is  of  size  128  x  128.  The  detail  coefficients  obtained  by  decomposing  the  image 
into  its  first  level  are  band-pass  filtered  by  a  filter  centered  at  48  cycles  per  image  with  cutoff 
frequencies  of  32  and  64  cpi  (octave).  See  figure  29. 


Cyci<»  par  imaga  Cydaa  par  imaga 

Figure  29.  Left:  Wavelet  high-pass  filter  centered  at  48  cpi,  Right:  Scale  low-pass  filter 
associated  with  the  hpf  on  left 

Applying  the  previous  steps  to  the  original  image  (defined  as  the  zeroth  level  approxima¬ 
tion  coefficient  set)  yielded  the  first  level  of  decomposition  coefficients.  To  obtain  additional 
levels  of  decomposition,  the  resulting  approximation  image  was  treated  as  the  image  and 
decomposed.  If  the  image  shown  in  figure  15  on  page  20  was  passed  through  the  algorithm 
outlined  in  figure  28,  each  decomposition  level  would  produce  an  approximation  image  related 
to  the  other  approximation  images  as  shown  in  figure  30. 

The  waveforms  used  for  reconstruction  (h  and  g)  were  calculated  from  their  decom¬ 
position  counterparts  (^  and  h)  using  equation  5;  that  is  g  was  used  to  generate  h  and  h  was 
used  to  generate  g.  These  new  waveforms  were  used  for  reconstruction  and  their  relationship 
is  shown  in  figure  3 1 . 

Reconstruction  of  the  original  image  was  the  mirror  operation  of  decomposition.  Fig¬ 
ure  32  shows  the  process  of  reconstructing  an  image  after  wavelet  decomposition.  The  exact 
steps  to  follow  for  reconstruction  of  an  image  decomposed  down  to  an  W  x  N  matrix  are: 
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Figure  30.  The  relationships  between  the  approximation  images  at  each  of  the  seven  levels 
of  decomposition  of  the  original  image  (largest  image  above).  Each  of  the  above 
images  would  serve  as  the  input  image  in  the  algorithm  of  figure  28.  These 
images  were  all  histogram  equalized  for  presentation  purposes. 
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Figure  31.  Reconstruction  scale  function  (h)  is  on  the  left,  reconstruction  wavelet  function 
(g)  is  in  the  middle,  while  their  relationship  is  shown  on  the  right. 
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Figure  32.  Reconstruction  of  an  image  after  a  single  level  of  decomposition 
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1 .  Upsample  the  rows  of  each  of  these  coefficient  sets  by  inserting  a  zeros  column  between 
every  column.  The  sizes  of  each  new  matrix  is  now  N  x  2N. 

2.  Convolve  the  rows  of  the  Approximation  and  Horizontal  coefficient  set  with  h  and  the 
rows  of  the  Vertical  and  Diagonal  coefficient  with  g. 

3.  Add  convolved  Approximation  image  and  Vertical  details  image  point-by-point.  Add 
convolved  Horizontal  details  image  and  Diagonal  details  image  point-by-point. 

4.  Upsample  the  columns  of  each  of  these  coefficient  sets  by  inserting  a  zeros  row  between 
every  row. 

5.  Generate  a  column  coefficient  set  by  convolving  the  columns  of  the  Approxima¬ 
tion/Vertical  information  image  with  g-,  generate  another  column  coefficient  set  by 
convolving  the  columns  of  the  Horizontal/Diagonal  information  image  with  h. 

6.  Add  these  two  remaining  matrices  point-by-point.  The  size  of  the  matrix  is  now 
2N  X  2N. 

The  resulting  image  after  the  above  operation  became  the  approximation  image  for  the  next 
higher  level  of  decomposition,  and  was  combined  in  an  identical  way  with  the  detail  coefficients 
at  the  next  higher  level.  For  example,  if  the  original  image  was  decomposed  down  to 
the  third  level,  reconstruction  would  first  involve  performing  figure  32  on  the  set  of  four 
matrices  representing  the  third  level  of  decomposition.  The  output  image  of  that  operation 
would  then  become  the  approximation  image  for  the  three  matrices  representing  the  second 
level  detail  coefficients  (horizontal,  vertical,  and  diagonal)  of  decomposition  and  figure  32 
would  be  performed  on  those  four  sets.  The  output  of  that  operation  would  become  the  first 
level  approximation  and  combined  with  the  first  level  details  to  produce  the  original  image. 
Appendix  C  shows  the  specific  Matlab  code  used  to  implement  this  biorthogonal  wavelet  set 
for  two  dimensional  images. 
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2.9  Headset 


The  headset  intended  to  be  the  output  for  image  enhancement  for  this  thesis  was  the 
PT-01;  a  display  system  which  could  be  raised  just  above  the  line  of  sight  and  not  obstruct 
peripheral  vision,  or  lowered  to  cover  the  eyes  completely.  It  contained  two  active  matrix 
420x230  pixel  screens.  Each  pixel  was  an  RGB  unit  containing  triads  of  three  cells  and 
occupying  less  than  3.2  arc  minutes  in  size.  Using  the  square  area  of  230x230  pixels  located 
in  the  center  of  the  displays,  the  total  area  would  correspond  to  12.3x12.3  degrees  of  arc. 
Therefore,  one  cycle  per  degree  equaled  12.3  cycles  per  image.  See  figure  33. 
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Figure  33.  Active  matrix  LCD  used  in  headset. 

2.10  Summary 

Millions  of  people  suffer  from  low  vision  due  to  macular  degeneration.  Knowledge  of 
a  low  vision  observers  visual  characteristics  could  perhaps  be  used  to  enhance  that  person’s 
view  of  the  world  so  that  face  recognition  could  be  improved.  Peli  had  good  success  with 
improving  the  recognition  abilities  of  low  vision  patients  using  adaptive  thresholding  and 
adaptive  filtering;  he  also  provided  a  simulation  which  could  give  those  with  normal  vision 


39 


an  insight  into  what  a  person  with  a  less  than  optimum  contrast  threshold  response  sees  when 
looking  at  a  given  image.  Lawton  and  Fine  demonstrated  that  compensation  filters  used  to 
enhance  words  improved  the  word  recognition  rate  in  low  vision  observers.  The  discrete 
wavelet  transform  (DWT)  was  seen  to  break  an  image  into  different  resolution  information 
while  preserving  location  information. 

The  next  chapter  outlines  some  new  wavelet  techniques  concerning  optimal  image 
enhancement  for  the  visually  impaired  and  a  new  wavelet  human  perception  simulation. 
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in.  Methodology 


The  research  area  outlined  here  was  new  to  AFIT.  This  chapter  presents  wavelet  analysis 
methods  that  will  be  compared  to  Peli’s  work  concerning  simulations  of  visual  perception  given 
a  particular  contrast  threshold  curve.  It  also  presents  image  enhancement  methods  that  use  a 
given  set  of  individual  contrast  perception  curves  to  influence  frequency  analysis  and  wavelet 
analysis.  A  new  method  using  suprathreshold  characteristics  of  a  patient  to  enhance  an  image 
is  also  presented  as  is  an  idea  concerning  localizing  the  enhancement. 

3.1  Wavelets 

In  a  practical  sense,  the  two-dimensional  Discrete  Wavelet  Transform  (DWT)  breaks 
an  image  into  a  set  of  images  corresponding  to  different  resolution  information  contained  in 
the  original  image.  For  example,  if  the  original  image  was  contained  in  a  128x128  array,  a 
wavelet  decomposition  would  generate  four  new  images  in  64x64  arrays,  four  in  32x32  arrays, 
four  in  16x16,  four  in  8x8,  four  in  4x4,  four  in  2x2,  and  four  in  1x1  arrays.  One  array  in  each 
set  of  four  would  correspond  to  low  frequency  information,  and  the  other  three  arrays  would 
correspond  to  high  frequency  information. 

The  method  by  which  an  image  was  decomposed  by  the  DWT  lent  itself  readily  to  the 
analysis  performed  on  contrast  in  complex  images  by  Peli  (19).  When  he  analyzed  his  image 
of  Cary  Grant  he  used  the  log  cosine  filters  shown  in  figure  16  on  page  21.  One  of  the  most 
useful  aspects  of  the  DWT  is  its  ability  to  locate  frequency  information  at  specific  locations 
in  the  original  image.  Peli’s  band-pass  filters  had  center  frequencies  at  powers  of  two  and 
were  one  octave  wide  (i.e.  the  filter  to  capture  the  information  around  16  cycles  per  image 
was  centered  at  16  cpi  with  cutoffs  at  12  and  24  cpi);  analogous  wavelet  filters  were  one 
octave  wide  with  cutoff  frequencies  at  powers  of  two  (for  example,  the  filter  to  capture  the 
information  around  12  cpi  was  centered  at  12  cpi  with  cutoffs  at  8  and  16  cpi).  Figures  34 
and  35  display  the  filters  of  the  parenthetic  examples  above.  It  is  important  to  remember  that 
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the  reason  for  the  upper  cutoff  of  the  high-pass  frequency  response  of  the  wavelet  is  that  the 
highest  frequency  is  limited  by  the  size  of  the  image  matrix. 


Figure  34.  Left;  High-pass  log-cosine  filter  centered  at  16  cpi,  Right:  Low -pass  log-cosine 
filter  associated  with  the  hpf  on  left. 


Figure  35.  Left:  High-pass  wavelet  filter  centered  at  12  cpi.  Right:  Low-pass  wavelet  filter 
associated  with  the  hpf  on  left. 


Since  the  discrete  wavelet  transform  used  circular  convolution  on  the  rows  and  columns 
during  decomposition  and  reconstruction,  edge  effects  became  apparent  whenever  the  detail 
coefficients  were  manipulated.  The  solution  used  in  this  thesis  was  to  implement  an  algorithm 
by  Cohen  and  do  the  DWT  on  a  “mirror  quad  image,”  such  as  the  one  shown  on  the  right  in 
figure  36,  instead  of  the  original  image  shown  on  the  left  (5).  This  operation  did  not  effect 
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Figure  36.  Mirror  quad  image  (size  256  x  256)  transformation  of  original  image.  The  image 
on  the  right  was  the  actual  input  into  the  algorithm  outlined  in  figure  28  on 
page  34. 


what  resolution  information  was  contained  in  each  of  the  coefficient  sets.  For  example,  for 
the  original  128  x  128  image,  the  first  approximation  yielded  a  low-pass  filtered  version  of 
the  original  with  cutoff  frequency  of  32  cycles  per  image.  Using  the  mirror  quad  operation, 
the  first  approximation  was  a  low-pass  filtered  version  of  the  mirror  quad  image  with  cutoff 
frequency  of  64  cycles  per  mirror  quad  image,  but  this  same  cutoff  was  still  32  cycles  per 
original  image.  The  decomposition  approximation  and  detail  coefficients  remained  in  the 
mirror  quad  format  during  any  analysis  of  the  information. 

An  important  consideration  after  reconstruction  was  that  the  reconstructed,  manipulated 
final  image  was  a  reconstructed  mirror  quad  image;  this  image,  though,  was  not  truly  a 
mirror  quad  image  due  to  the  various  manipulations  encountered  between  decomposition  and 
reconstruction.  To  maintain  manipulated  information  from  each  quadrant  of  the  final  image, 
the  displayed  image  was  an  average  of  the  quadrants  of  the  reconstructed  mirror  quad  image. 
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3.2  Wavelet  Simulations 


This  section  outlines  how  the  bi-orthogonal  wavelet  set  was  used  to  parallel  the  work 
done  by  Peli  concerning  simulation  of  image  perception  by  a  subject  with  a  given  contrast 
sensitivity  curve  (19). 

The  first  step  in  the  process  was  to  generate  a  set  of  images  filtered  from  the  original 
using  octave-wide  band-pass  filters  and  corresponding  low-pass  filters  which  captured  the 
information  under  the  band-pass  filters.  Wavelet  decomposition  levels  provided  the  ideal 
relationship  in  that  there  were  coefficient  sets  containing  detail  information  about  the  image 
at  different  resolutions  which  could  be  compared  to  the  approximation  coefficient  set  at  that 
level.  For  example,  the  third  level  details  (shown  in  figure  37)  were  each  point-by-point 


Figure  37.  The  first  image  above  contains  the  third  level  vertical  detail  coefficients,  the 
second  is  the  horizontal  detail  coefficient  set,  and  the  last  image  is  the  diagonal 
detail  coefficient  set. 


Figure  38.  The  image  above  is  the  third  level  decomposition  approximation  coefficient  set. 

divided  by  the  approximation  coefficients  in  figure  38,  the  end  result  being  three  contrast 
images  representing  vertical,  horizontal,  and  diagonal  contrasts  at  those  pixel  positions  at  that 
level  of  decomposition.  This  first  step  required  some  subtle  manipulations  while  decomposing 
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the  original  image.  The  third  level  decomposition  sets  would  appear  as  in  figures  37  and  38 
(discounting  the  mirror  quad  operation);  however,  it  was  inappropriate  to  use  these  for  the 
simulation  due. to  the  shifting  and  downsampling  imposed  by  the  wavelet  decomposition. 
Referring  back  to  figure  27,  notice  that  the  wavelet  function  (p)  is  not  symmetric  about  the  same 
point  as  the  scaling  function  (h).  Whenever  an  image  is  being  decomposed,  the  information 
will  not  be  shifted  by  the  convolution  with  the  scale  function  since  it  is  symmetric  about  zero; 
however,  the  convolution  with  the  wavelet  function  produces  a  shift  of  the  information  since 
the  wavelet  function  is  symmetric  about  the  first  index  to  the  left  of  zero.  In  order  to  parallel 
Dr.  Peli’s  log-cosine  filter  analysis,  each  of  the  detail  coefficient  sets  at  each  level  were  to 
be  divided  by  the  approximation  image  at  that  level;  this  division  would  yield  a  “resolution 
limited  contrast  image.”  Since  the  convolution  with  the  wavelet  function  would  shift  the 
resulting  detail  coefficients  with  respect  to  the  rows  and  columns  convolved  with  the  scale 
function,  this  shifting  aspect  had  to  be  eliminated  during  the  resolution  limited  contrast  image 
generation  stage.  In  order  to  eliminate  this  shifting  during  convolution,  the  wavelet  was  forced 
to  center  at  zero.  Figure  39  shows  the  new,  modified  scale  and  wavelet  vectors  used  in  the 
wavelet  simulation  decomposition  and  figure  40  shows  how  the  un-shifted  detail  coefficient 
sets  were  generated.  Compare  figure  39  to  figure  27  on  page  33  and  notice  that  the  functions 


Figure  39.  Scale  function  (at  left)  and  un-shifted  wavelet  function  (middle)  used  in  wavelet 
Simulation.  Their  relationship  is  shown  in  the  plot  on  the  right.) 

are  the  same,  only  the  indexing  is  different.  This  new  wavelet  decomposition  set  was  then 
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convolved  with  the  rows  and  columns;  this  allowed  the  resolution  limited  contrast  images  to 
be  calculated  with  the  positions  of  all  the  decomposition  images  correctly  aligned  relative  to 
one  another. 
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Figure  40.  Generation  of  un-shifted  detail  coefficients. 


The  result  of  using  the  algorithm  of  figure  40  was  that  the  first  level  of  decomposition 
yielded  four  matrices  that  were  the  same  size  as  the  original  image  matrix,  and  those  four 
matrices  were  used  to  compute  the  resolution  limited  contrast  images. 

These  resolution  limited  contrast  images  were  compared  point-by-point  to  the  proper 
threshold  value  (interpolated  from  the  contrast  threshold  curve  in  figure  13  on  page  18  used 
for  the  other  simulations).  If  the  value  was  at  or  above  threshold,  the  detail  coefficient  at  that 
location  was  left  alone.  Otherwise,  if  the  value  was  below  threshold,  the  detail  coefficient 
at  that  location  was  set  to  zero.  For  example,  a  set  of  resolution  limited  contrast  images 
corresponding  to  12  cycles  per  image  was  generated,  using  the  coefficient  sets  in  the  third 
level  decomposition  (shown  in  figure  37  and  38)  by  dividing  each  of  the  three  third  level  detail 
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coefficient  sets  point-by-point  by  the  third  level  approximation  image.  These  three  contrast 
images  were  point-by-point  compared  to  the  value  interpolated  at  12  cpi  from  the  plot  for  the 
person  with  20/136  vision  shown  in  figure  13  on  page  18  (that  threshold  value  was  .4096). 
If  the  value  at  that  coordinate  position  in  the  detail  contrast  image  was  below  threshold,  the 
pixel  value  in  the  corresponding  detail  coefficient  matrix  was  set  to  zero.  This  represented  a 
loss  of  detail  at  that  coordinate  position.  If  the  value  at  that  coordinate  position  in  the  detail 
contrast  image  was  above  threshold,  the  value  in  the  corresponding  detail  coefficient  set  was 
left  untouched. 

It  was  at  this  point  in  the  process  where  downsampling  occurred.  The  critical  shift  that 
was  eliminated  previously  now  had  to  be  reinserted  because  it  was  vital  to  the  reconstruction 
process.  This  algorithm  did  the  shifting  during  the  downsampling  process  instead  of  during 
the  convolution  with  the  wavelet  function;  the  downsampling  of  any  column  or  row  convolved 
with  the  scaling  function  took  the  even  indices  starting  with  index  two,  while  the  downsampling 
of  any  column  or  row  convolved  with  the  wavelet  function  took  the  odd  indices  starting  with 
index  three  (the  last  index  being  index  one  due  to  circular  convolution).  Figure  41  shows 
how  the  downsampling  was  used  to  implement  the  shifting.  This  preserved  the  information 
storage  aspect  of  the  wavelet  decomposition  and  allowed  reconstruction  to  occur  in  the  normal 
fashion  as  shown  in  figure  32  on  page  37. 

This  process  was  repeated  at  all  resolution  levels  down  to  the  smallest  resolution  level 
with  detail  coefficients  which  had  maximum  resolution  of  4  cycles  per  image  (which  resulted 
in  a  4  X  4  image  after  downsampling),  and  then  reconstructed  in  the  standard  fashion  of  the 
biorthogonal  wavelet,  as  shown  in  figure  32  on  page  37.  Appendix  D  contains  the  Matlab 
code  used  for  these  simulations. 

3.3  Wavelet  Enhancement 

In  order  to  maintain  the  flow  of  the  work  already  accomplished  by  Dr.  Peli,  the  first 
wavelet  enhancement  technique  explored  here  paralleled  his  modified  pre-emphasis  filter 
enhancement  technique.  The  pre-emphasis  model  boosted  the  highest  frequencies  the  subject 
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Figure  41.  Method  by  which  the  downsampling  was  used  to  achieve  the  shifting  normally 
accomplished  by  the  wavelet. 

could  see  by  the  inverse  of  the  normalized  threshold  measured  at  that  frequency  for  that  person. 
For  the  wavelet  case,  the  boosting  was  done  to  the  highest  detail  coefficients  the  person  would 
be  able  to  discern.  In  Dr.  Peli’s  paper  in  which  the  pre-emphasis  model  is  presented  (24),  the 
example  centered  the  boosting  frequency  at  16  cycles  per  image.  The  wavelet  model  presented 
here  boosted  the  frequencies  centered  at  12  cycles  per  image;  this  new  center  frequency  was 
determined  by  the  nature  of  the  discrete  wavelet  decomposition.  Figures  37  and  38  show  the 
third  level  detail  coefficients  and  the  third  level  approximation  coefficients.  With  the  detail 
coefficients  generated  down  to  the  desired  level  to  be  boosted  (here  that  was  level  three),  each 
detail  coefficient  set  was  multiplied  by  whatever  the  enhancement  algorithm  called  for  at  each 
specific  resolution  based  on  an  individuals  perception  characteristics  (here  the  third  level  detail 
coefficients  were  multiplied  by  a  factor  of  five).  Each  of  the  three  detail  levels  (horizontal, 
vertical,  and  diagonal)  were  multiplied  by  the  same  number,  since  the  combination  of  the  three 
corresponded  to  a  band-pass  filtered  image. 
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3.4  Using  Individual  Suprathreshold  Contrast  Characteristics  of  a  Patient  for  Enhancement 

Peli’s  early  work  in  digital  image  enhancement  relied  heavily  on  boosting  frequencies 
most  useful  for -face  recognition  without  regard  for  a  particular  person’s  visual  characteristics 
(26,  22,  24).  His  and  other’s  later  experimentation  on  tailoring  enhancement  to  the  individual 
was  based  on  measuring  contrast  threshold  characteristics  of  a  patient  and  using  that  infor¬ 
mation  for  enhancement  (15,  19,  21,  24,  7).  Rubin  argued  (based  on  testing  of  both  young 
and  old  observers)  that  the  contrast  sensitivity  of  a  low  vision  observer  may  not  in  fact  predict 
recognition  ability  (31);  here  it  is  proposed  that  suprathreshold  characteristics  of  a  subject 
would  be  more  useful  for  image  enhancement.  The  contrast  in  facial  images  is  usually  above 
the  threshold  of  normal  observers.  Using  contrast  sensitivity,  which  is  based  on  contrast 
threshold,  for  enhancement  may  not  optimize  the  information  in  the  image  that  the  subject 
can  see.  Suprathreshold  contrast  perception  curves,  close  to  threshold,  tend  to  be  shaped 
the  same  as  the  threshold  curve;  but  as  the  contrast  increases,  the  relationship  between  the 
perceived  suprathreshold  contrasts  changes.  If  the  suprathreshold  contrast  perception  of  the 
subject  was  used  to  enhance  the  image  instead  of  the  threshold  contrast  perception,  this  would 
more  efficiently  enhance  those  frequencies  dominant  in  the  image.  The  remainder  of  this 
section  describes  a  new  method  of  image  enhancement  tailored  to  an  individual  subject  using 
measured  normalized  suprathreshold  curves  such  as  those  shown  in  figure  14  on  page  19  (2). 

First,  the  contrast  information  in  the  image  had  to  be  determined.  Since  Cannon's 
suprathreshold  matching  experiments  (see  figure  14  on  page  19)  indicated  that  4  cycles 
per  degree  was  the  most  sensitive  frequency  of  those  normal  observers  measured  during 
that  particular  study,  that  frequency  was  the  “reference  frequency.”  Different  measurement 
techniques  may  yield  different  suprathreshold  contrast  curves,  thereby  determining  the  actual 
reference  frequency.  The  technique  used  to  determine  the  contrast  at  that  frequency  took 
advantage  of  the  band-limited  contrast  concept  outlined  in  his  paper  on  contrast  in  complex 
images  (19).  The  image  was  band-pass  filtered  with  Peli’s  log-cosine  filter  corresponding  to 
4  cpd  to  create  the  bpf  image;  low-pass  filtered  with  the  low-pass  filter  to  complement  that 
band-pass  filter  (analogous  to  the  filters  shown  in  figure  16  on  page  21),  and  the  former  image 
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point-by-point  divided  by  the  latter.  This  operation  created  a  reference  contrast  image  for  the 
frequency  most  detectable  to  the  normal  observer.  The  global  reference  contrast  of  the  image 
was  then  deterrnined  to  be  the  root  mean  square  of  the  pixel  values  in  the  image.  This  global 
reference  contrast  of  the  4  cpd  information  in  the  image  determined  which  of  the  suprathreshold 
curves  to  use.  For  example,  if  figure  14  on  page  19  was  the  set  of  curves  generated  by  testing 
a  low-vision  subject  (even  though  these  particular  curves  were  generated  using  the  average  of 
a  number  of  normal  observers),  and  the  rms  value  of  the  reference  contrast  image  fell  closest 
to  .01  (as  opposed  to  .003  or  .03),  the  curve  chosen  for  enhancement  purposes  would  be  the 
first  curve  up  from  the  threshold  curve  (where  the  contrast  of  4  cpd  was  .01). 

Now  that  the  proper  curve  was  chosen,  the  frequency  manipulation  to  be  performed  was 
determined  by  the  suprathreshold  curve  generated  by  normalizing  the  curve  for  the  subject  to 
a  curve  generated  by  a  number  of  normal  observers.  This  new  curve  was  the  “Normalized 
Contrast  Suprathreshold  Curve”  (NCSC).  Its  inverse  was  used  as  an  inverse  filter,  much  in  the 
same  manner  as  what  Lawton  did  in  her  experiments  concerning  improving  reading  rates  in 
low-vision  observers  (15). 

3.5  Localization 

Unless  an  observer  is  trying  to  recognize  the  face  of  a  person  very  distant  from  him,  the 
fovea  of  the  observer  will  not  usually  encompass  the  entire  visage  of  the  subject  of  interest. 
The  fovea  covers  approximately  one  to  two  degrees  of  arc  in  the  visual  field  of  view  (1),  and 
for  this  reason  people  “saccade”  and  scan  the  face  of  interest  to  get  the  “whole  picture.”  This 
implies  that,  at  any  given  time  during  this  process,  the  observer  is  using  his  fovea  on  only  a 
small  part  of  the  facial  image;  the  rest  of  the  face  is  still  in  view,  but  being  perceived  by  the 
less  visually  acute  retinal  periphery.  These  facts  lead  to  the  problem  of  how  best  to  locally 
enhance  a  face  for  individualized  perception. 

Global  manipulation  of  an  image  wastes  many  precious  resources:  the  finite  dynamic 
range  of  the  screen  limits  the  possible  gain  of  high  frequencies  (21);  enhancement  algorithms 
work  more  slowly  when  processing  the  whole  scene;  and  the  person  is  already  placing  an 
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object  of  interest  in  their  area  of  highest  visual  acuity  (24).  Since  only  a  certain  amount  of 
information  is  being  processed  at  a  given  time  as  an  observer  peers  at  a  facial  image,  one 
method  of  optimization  could  be  to  break  the  image  into  “chunks”  and  enhance  each  chunk 
separately  and  simultaneously.  The  weakness  with  this  method  is  that  the  observer  may  not 
want  to  look  only  at  predetermined  discrete  chunks  of  an  image,  but  may  want  to  choose 
which  chunks  to  focus  upon,  saccading  the  image  and  gathering  information  at  each  stopping 
point.  This  fact  implies  that  the  optimum  solution  would  be  to  enhance  the  image  only  in  the 
area  of  interest  to  the  observer. 

Retina  motion  trackers  have  been  developed  (11,  38)  which  can  tell  precisely  where 
each  eye  is  looking  at  any  given  instant.  Since  persons  with  macular  defects  tend  to  place  an 
object  of  interest  just  outside  the  dead  area  of  the  fovea  (24),  if  that  low  vision  subject  were 
instructed  to  look  at  a  point  of  interest  corresponding  to  one  or  two  degrees  of  arc  so  as  to 
see  that  point  in  the  most  optimal  fashion,  the  motion  tracker  could  be  aligned  to  this  new 
“pseudo-retina.”  The  enhancement  of  a  given  image  would  then  be  tied  directly  to  what  the 
subject  wanted  enhanced  (i.e.  where  the  subject  was  looking). 

3.6  Metric  of  Enhancement  Ejfect 

The  intrinsic  complexity  and  adaptability  of  the  human  visual  system  has  caused  more 
than  a  little  consternation  in  the  image  processing  community.  Image  quality  is  a  subjective 
term  that  refuses  to  be  tied  down  to  a  single  set  of  numbers  or  equations  which  may  allow 
a  solid  metric  by  which  image  manipulation  results  can  be  evaluated.  Images  of  the  natural 
world  contain  many  frequencies  and  intensities  for  the  human  eye  to  process;  add  to  this  the 
restricted  field  of  view  imposed  by  a  headset  and  the  limited  dynamic  range  of  the  output 
screens  of  the  image  enhancement  system  being  viewed  by  the  subject,  and  the  problem 
becomes  even  more  difficult. 

Here,  one  metric  of  enhancement  effect  is  proposed.  Although  image  enhancement 
usefulness  for  any  given  person  can  only  be  determined  experimentally,  the  enhancement 
effect  metric  proposed  here  was  based  on  the  known  edge  detection  properties  incorporated  in 
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the  average  human  visual  system.  First,  the  assumption  was  made  that  a  simulation  of  human 
visual  perception  based  on  an  individual’s  contrast  sensitivity  is  representative  of  what  that 
individual  perceives  when  looking  at  an  image.  Both  the  Peli  and  wavelet  simulations  were 
used.  Second,  each  enhancement  to  be  evaluated  was  run  through  that  simulation.  Third, 
an  edge  detection  algorithm  was  applied  to  the  original  image  and  the  output  images  of  the 
simulation.  Finally,  the  edge  detected  enhanced  images  were  compared  to  the  edge  detected 
original  image  by  cross  correlating  the  two  edge-detected  images;  the  desired  result  being  an 
image  perceived  by  the  individual  that  is  perceived  to  be  as  close  as  possible  to  the  original. 
The  edge  detection  convolution  mask  used  was  the  Laplacian  mask  (9).  This  mask,  shown  in 
figure  42,  was  two-dimensionally  convolved  with  the  image  of  interest.  This  resulting  image 


Figure  42.  Laplacian  convolution  mask  used  for  edge  detection. 


was  a  measure  of  the  effective  edges  in  the  image. 

As  an  example,  using  this  metric  to  measure  the  effective  edge  enhancement  of  the 
pre-emphasis  filter  enhancement  outlined  in  section  2.7.5,  the  first  step  involved  detecting 
the  edges  in  the  original  image  shown  on  the  left  in  figure  43.  The  resulting  edge-detected 
image  is  shown  on  the  right.  The  image  on  the  left  in  figure  44  is  the  pre-emphasis  filter 
enhanced  image  after  being  run  through  the  Peli  simulation,  while  the  image  on  the  right  is 
the  edge-detected  version  of  the  image  on  the  left.  The  similarity  measure  was  computed 
by  performing  point-by-point  multiplication  of  the  edge-detected  original  (on  the  right  of 
figure  43)  and  the  edge-detected  enhanced-perceived  image  (on  the  right  of  figure  44,  then 
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Figure  43.  Left:  original  image.  Right:  edge-detected  original  image. 


Figure  44.  Left:  edge-detected  pre-emphasis  filtered  after-simulation  image.  Right:  edge- 
detection  resultant  of  image  on  left. 

summing  the  resulting  matrix  values  to  yield  a  scalar.  This  scalar  was  divided  by  the  scalar 
obtained  from  comparing  the  original  to  itself  (autocorrelation).  This  operation  was  intended 
to  normalize  the  data  between  zero  and  one.  If  zero,  the  resulting  edge  analysis  yielded  no 
similarity  to  the  original  image’s  edge  analysis;  if  one,  the  edges  were  identical  to  the  original. 
Our  example  above  yielded  a  correlation  value  of  0.1885. 

3.7  Summary 

A  biorthogonal  wavelet  set  was  used  to  simulate  what  a  person  with  a  given  contrast 
sensitivity  function  sees  when  looking  at  an  image.  This  same  wavelet  set  was  used  to  enhance 
an  image  in  a  manner  analogous  to  Peli’s  modified  pre-emphasis  filter  enhancement.  A  new 
method  for  image  enhancement  based  on  a  subject’s  measured  suprathreshold  characteristics 
was  outlined.  Finally,  the  concept  was  proposed  that  the  image  should  be  enhanced  only  in 
the  area  of  interest  to  the  observer.  The  results  of  the  wavelet  simulation  and  enhancement  are 
presented  and  analyzed  in  the  next  chapter. 


IV.  Results 


This  chapter  outlines  the  results  of  trying  the  various  methods  discussed  in  the  previous 
chapter. 

4.1  Simulation  of  Macular  Degeneration  Using  Wavelets 

When  using  the  wavelet  scheme  to  simulate  the  vision  of  the  person  with  the  threshold 
curve  corresponding  to  acuity  of  20/136,  the  image  on  the  left  in  figure  45  was  operated  on 
and  the  resulting  image  was  as  shown  in  the  middle.  The  image  in  the  middle  is  analogous 


Figure  45.  Left:  Original  image.  Middle:  Wavelet  simulation  of  what  person  with  threshold 
curve  corresponding  to  20/136  visual  acuity  shown  on  page  1 8  sees  when  looking 
at  image  on  left.  Right:  Peli  simulation  of  what  that  same  person  sees. 

to  the  image  created  by  the  Peli  simulation  shown  on  the  right.  Notice  that  the  wavelet 
simulation  indicates  that  less  information  is  lost  than  is  indicated  by  the  Peli  simulation. 
This  is  due  to  the  fact  that  the  reconstruction  process  performed  by  the  wavelet  analysis  is 
different  than  the  reconstruction  process  performed  by  Peli’s  analysis.  In  Peli’s  simulation, 
each  image  corresponding  to  a  band-pass  filtered  original  image  was  manipulated;  these 
resulting  images  were  then  simply  added  together  to  form  the  image  on  the  right.  In  the 
wavelet  simulation,  each  image  corresponding  to  a  detail  coefficient  set  (horizontal,  vertical, 
or  diagonal)  was  manipulated;  these  resulting  coefficient  sets  were  then  recombined  at  each 
level  in  the  standard  discrete  wavelet  reconstruction  fashion;  that  is  to  say,  many  convolutions 
were  performed  during  the  reconstruction  phase  that  affected  the  final  output  image. 
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When  the  curve  representing  normal  vision  in  figure  13  was  used  with  the  previous 
wavelet  simulation  on  the  presented  image  on  the  left  in  figure  46,  the  image  shown  in  the 
middle  is  the  result.  Notice  that  even  less  distortion  is  evident  than  the  same  Peli  simulation 


Figure  46.  Left:  Original  image.  Middle:  Wavelet  simulation  of  what  person  with  threshold 
curve  representing  normal  vision  shown  on  page  18  sees  when  looking  at  image 
on  left.  Right:  Peli  simulation  of  what  that  same  person  sees. 


shown  on  the  right. 

Whenever  the  image  enhanced  using  adaptive  filtering  shown  on  the  left  in  figure  47  is 
passed  through  the  wavelet  simulation  discussed  above,  the  result  is  as  shown  in  the  middle 
image.  The  wavelet  simulation  indicates  that  the  adaptive  filtering  enhancement  boosts  the 


Figure  47.  Left:  Image  enhanced  using  adaptive  filtering.  Middle:  Wavelet  simulation  of 
what  person  with  threshold  curve  corresponding  to  20/136  visual  acuity  shown 
on  page  18  sees  when  looking  at  image  on  left.  Right:  Peli  simulation  of  what 
that  same  person  sees. 


information  available  to  the  low-vision  observer.  The  image  in  the  middle  is  analogous  to  the 
image  presented  on  the  right.  The  wavelet  simulation  seems  to  once  again  indicate  that  less 
information  is  lost  than  the  Peli  simulation  would  indicate;  this  being  due  to  the  reconstruction 
process  in  the  wavelet  analysis. 
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The  result  of  running  the  image  enhanced  using  Lawton’s  compensation  filter  from 
section  2.7.4  (shown  at  left  in  figure  48)  through  the  wavelet  simulation  is  shown  in  the  middle 
image.  Notice  once  again  that  the  middle  image  clearly  contains  more  information  than  the 


Figure  48.  Left:  Image  enhanced  using  compensation  filtering.  Middle:  Wavelet  simulation 
of  what  person  with  threshold  curve  corresponding  to  20/136  visual  acuity  shown 
on  page  18  sees  when  looking  at  image  on  left.  Right:  Peli  simulation  of  what 
that  same  person  sees. 

Peli  simulation  image  on  the  right  indicates. 

When  the  image  enhanced  using  the  modified  pre-emphasis  filter  (shown  on  the  left  in 
figure  49)  was  passed  through  the  wavelet  simulation  of  what  a  person  with  20/136  vision 
might  see  when  looking  at  that  image,  the  result  was  as  shown  in  the  middle  image.  This 


Figure  49.  Left:  Image  enhanced  using  modified  pre-emphasis  filtering,  Middle:  Wavelet 
simulation  of  what  person  with  threshold  curve  corresponding  to  20/136  visual 
acuity  shown  on  page  18  sees  when  looking  at  image  on  left.  Right:  Peli 
simulation  of  what  that  same  person  sees. 


image  is  analogous  to  the  image  shown  on  the  right.  Here,  the  simulated  image  in  the  middle 
is  virtually  identical  to  the  enhanced  image  on  the  left  due  to  the  reconstruction  process. 
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4.2  Wavelet  Enhancement 


When  the  third  level  of  details  were  multiplied  by  a  factor  of  five,  and  the  image 
reconstructed,  the  resulting  image  was  as  (shown  on  the  left  in  figure  50),  with  the  image  in 
the  middle  in  that  same  figure  representing  what  a  person  with  20/136  visual  acuity  sees  when 
viewing  the  image  on  the  left  according  to  the  wavelet  simulation.  The  image  on  the  left  is 


Figure  50.  Left:  Original  image  enhanced  by  multiplying  the  third  level  detail  coefficients 
by  five.  Middle:  Wavelet  simulation  of  what  person  with  threshold  curve  corre¬ 
sponding  to  20/136  visual  acuity  shown  on  page  18  sees  when  looking  at  image 
on  left.  Right:  Peli  simulation  of  what  that  same  person  sees. 

analogous  to  the  modified  pre-emphasis  filter  shown  in  figure  26  on  page  30.  Comparing  the 
images  on  the  right  of  these  two  figures,  the  pre-emphasis  filter  seems  to  preserve  some  of  the 
fine  details  associated  with  the  original  image  very  well,  while  the  wavelet  enhanced  version 
preserves  different  information  about  the  details. 

Diverging  slightly  from  the  specific  incarnation  of  Peli’s  modified  pre-emphasis  filter 
outlined  in  his  paper  on  image  enhancement  for  the  visually  impaired  (24)  by  multiplying  the 
third  level  coefficients  by  two  instead  of  five  results  in  the  image  on  the  left  in  figure  5 1 .  This 
enhancement  proved  more  subtle  as  can  be  seen  in  the  image  in  the  middle.  The  sharp  edges 
present  in  figure  50  have  been  reduced  at  the  possible  expense  of  detail. 

4.3  Metric  of  Enhancement  Effect 

The  use  of  edge  detection  analysis  on  enhancement  effectiveness  using  the  Peli  sim¬ 
ulations  resulted  in  values  as  shown  in  table  3,  while  the  use  of  edge  detection  analysis  on 
enhancement  effectiveness  using  the  wavelet  simulations  resulted  in  values  as  shown  in  ta- 
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Figure  5 1 .  Left:  Original  image  enhanced  by  multiplying  the  third  level  detail  coefficients 
by  two.  Middle:  Wavelet  simulation  of  what  person  with  threshold  curve  corre¬ 
sponding  to  20/136  visual  acuity  shown  on  page  18  sees  when  looking  at  image 
on  left.  Right:  Peli  simulation  of  what  that  same  person  sees. 

ble  4.  Both  analyses  clearly  indicate  that  any  of  the  outlined  manipulations  of  the  image  boost 
the  edge  energy  being  perceived  by  the  subject  beyond  what  was  obtained  with  no  enhance¬ 
ment.  The  values  calculated  that  were  above  one  occurred  with  thresholding  and  adaptive 
filtering;  this  may  indicate  too  much  edge  enhancement  since  subjects  have  complained  about 
the  thresholding  technique  making  the  image  look  too  cartoonish.  Note  also  that  the  wavelet 
enhancement  results  were  roughly  the  same  as  the  modified  pre-emphasis  filter  enhancement 
results.  This  modified  pre-emphasis  filter  is  currently  what  is  being  used  in  the  LVES  at 
Johns-Hopkins. 


Enhancement 

Similarity 

Original  image 

1 

No  enhancement 

0.0092 

Thresholding 

0.1019 

Adaptive  filter 

0.2003 

Compensation  filter  (no  histogram  equalization) 

0.0411 

Compensation  filter  (histogram  equalized) 

0.1953 

Pre-emphasis  filter 

0.1885 

Wavelet  enhance(gain=2) 

0.0660 

Wavelet  enhance(gain=3) 

0.0949 

Wavelet  enhance(gain=4) 

0.0972 

Wavelet  enhance(gain=5) 

0.1106 

Table  3.  Applying  edge  detection  metric  to  perceived  (Peli  simulation)  enhanced  images. 
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Enhancement 

Similarity 

Original  image 

1 

No  enhancement 

0.3154 

Thresholding 

1.4091 

Adaptive  filter 

1.3047 

Compensation  filter  (no  histogram  equalization) 

0.4633 

Compensation  filter  (histogram  equalized) 

0.3378 

Pre-emphasis  filter 

0.3464 

Wavelet  enhance(gain=2) 

0.4183 

Wavelet  enhance(gain=3) 

0.4464 

Wavelet  enhance(gain=4) 

0.4340 

Wavelet  enhance(gain=5) 

0.4237 

Table  4.  Applying  edge  detection  metric  to  perceived  (wavelet  simulation)  enhanced  images. 
4.4  Summary 

The  result  of  using  a  biorthogonal  wavelet  set  to  simulate  the  perception  of  a  person 
with  a  given  set  of  visual  perception  characteristic  curves  indicates  that  the  details  of  an 
image  may  not  be  as  suppressed  as  the  simulation  proposed  by  Peli  indicates.  The  new 
enhancement  algorithms  using  wavelets  are  seen  to  be  comparable  to  the  enhancements 
previously  explored  by  Peli  and  others.  The  enhancements  occur  primarily  along  the  horizontal 
and  vertical  directions,  which  is  reasonable  due  to  the  nature  of  the  wavelet  decomposition 
and  reconstruction.  Figure  52  displays  many  of  the  enhancement  techniques  and  the  effects 
of  running  the  images  through  the  simulations.  The  next  chapter  discusses  the  implications  of 
the  wavelet  simulations  and  enhancements. 


59 


Figure  52.  Summary  of  various  enhancements.  The  left  column  contains  enhancements,  the 
middle  column  contains  the  image  at  the  left  run  through  the  Peli  simulation,  the 
right  column  contains  the  image  at  the  far  left  run  through  the  wavelet  simulation 
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V.  Conclusions 


5.1  Analysis  of  Algorithms  Reproduced 

The  previous  work  in  the  area  of  image  interpretation  and  enhancement  for  the  visually 
impaired  concentrated  primarily  on  the  manipulation  of  the  Fourier  spectrum.  One  desire 
of  past  research  has  been  to  use  the  visual  characteristics  of  the  observer  as  a  guide  for 
enhancement,  since  those  visual  characteristics  should  be  able  to  be  compared  to  a  normal 
observer’s  visual  characteristics  and  the  image  enhanced  somehow  based  on  this  knowledge. 
It  has  been  conjectured  that  the  visual  system  of  a  subject  could  be  understood,  at  least  in 
part,  by  non-invasively  measuring  the  contrast  threshold  and  suprathreshold  characteristics  of 
that  subject.  This  research  reproduced  Peli’s  perception  simulation  model  using  the  contrast 
threshold  response  of  an  observer  to  simulate  the  appearance  of  an  image  to  that  low-vision 
observer.  This  simulation  was  based  on  filtering  characteristics  of  the  human  visual  system 
and  its  inability  to  perceive  those  frequency  components  below  threshold  contrast.  Various 
enhancement  schemes  from  the  literature  were  reproduced,  most  of  which  manipulated  the 
Fourier  spectrum  of  the  image  to  make  it  more  recognizable  to  the  low-vision  viewer.  The 
compensation  filter  used  by  Lawton  to  enhance  letters  (15)  was  used  to  enhance  faces.  This 
research  demonstrated  it  was  a  feasible  tool  for  image  enhancement.  Peli’s  modified  pre¬ 
emphasis  filter  enhancement,  which  is  being  used  in  the  LVES,  was  reproduced  and  seen  to 
also  enhance  the  details  significantly. 

A  practical  application  of  the  discrete  wavelet  transform  was  implemented  and  dis¬ 
cussed,  as  was  its  relationship  to  Fourier  analysis.  'Wavelet  analysis  was  seen  to  break  an 
image  into  different  global  resolution  information  sets  which  could  be  manipulated  and  re¬ 
constructed. 

5.2  Analysis  of  Original  Algorithms 

Here  was  proposed  a  new  simulation  of  human  perception  based  on  wavelet  analysis. 
This  simulation  using  wavelets  of  what  a  person  with  a  given  threshold  curve  sees  when 
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looking  at  an  image  yielded  different  results  than  the  analogous  simulation  using  the  log 
cosine  filters  proposed  by  Peli.  Most  noticeable  was  the  fact  that  more  details  were  present 
in  the  resulting  wavelet  simulation  images  than  in  the  Peli  simulation  images.  Using  wavelet 
analysis  to  simulate  human  visual  perception  provided  a  different  choice  for  a  basis  set,  one 
that  is  supported  by  recent  work  concerning  mammalian  cortical  processing  (12,  13,  35). 
The  overriding  concern  with  all  of  these  simulations  is  that  human  perception  is  still  poorly 
understood.  If  wavelet  analysis  truly  does  model  human  visual  perception,  then  the  next 
logical  step  is  to  determine  which  set  of  coefficients  the  brain  uses  to  process  that  information. 
The  wavelet  simulation  implies  that  the  various  enhancement  algorithms  presented  here  are, 
in  fact,  perceived  by  the  subject.  One  test  to  measure  the  value  of  the  wavelet  simulation 
would  be  to  have  a  person  with  normal  vision  in  one  eye  and  low  vision  in  the  other  view  an 
image  one  eye  at  a  time  and  give  an  opinion. 

Wavelets  were  used  to  provide  a  new  image  enhancement  procedure  by  manipulating 
the  decomposition  detail  coefficients.  The  new  metric  to  analyze  the  enhancements  based 
on  edge  similarity  measures  showed  that,  at  least  in  the  regard  of  edge  enhancement,  the 
perceived  wavelet  enhancements  were  matching  the  original  image  at  about  the  same  level 
as  the  perceived  modified  pre-emphasis  filter  enhancements  proposed  by  Peli.  The  true 
usefulness  of  the  wavelet  coefficients  as  enhancement  tools  can  only  be  determined  by  the 
testing  process.  The  only  sure  statement  that  can  be  made  about  image  enhancement  for  the 
visually  impaired,  without  first  having  a  better  understanding  of  the  human  visual  system, 
is  that  individual  subject  preferences  will  be  the  deciding  factor  when  choosing  the  optimal 
enhancement  technique.  Wavelet  enhancements  would  provide  different  information  to  the 
low-vision  observer  and  could  prove  to  be  superior  to  other  enhancement  techniques  for  some 
patients. 

5.3  Recommendations  for  Follow-On  Research 

The  most  obvious  area  that  was  not  developed  in  this  research  was  actual  testing  of  the 
wavelet  enhancement  techniques.  The  steps  necessary  to  do  this  would  include: 
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•  Measure  the  output  of  the  headset  with  a  photometer  and  develop  a  curve  so  that  the 
non-linearity  of  the  screens  can  be  accounted  for. 

•  Develop  a  digital  to  analog  conversion  scheme  so  that  the  output  displays  of  the 
headset  do  not  display  Mach  bands  at  the  very  low  (threshold)  contrast. 

•  Test  visually  impaired  subjects  for  contrast  sensitivity  and  contrast  suprathreshold 
characteristics. 

•  Incorporate  the  suprathreshold  contrast  enhancement  scheme  outlined  here  with  both 
Fourier  techniques  and  wavelet  techniques. 

•  Develop  software  to  take  advantage  of  the  oculometer  which  tracks  retinal  movement 
so  that  just  the  area  of  the  retina  in  the  subject  receives  an  enhanced  image  segment. 

•  Test  wavelet  enhancement  schemes  using  bi-orthogonal  wavelets  as  well  as  other 
symmetric  wavelets. 
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Appendix  A.  Peli  Simulation 


The  simulation  of  what  a  person  with  a  give  contrast  threshold  sees  was  accomplished 
with  the  m-file  scotoma.m: 


%This  function  will  simulate  what  a  person  with  central  scotoma  sees  (outputtmage)  when  looking  at  a  picture  (inputimage) . 

function  outputimage-scotoma(inputimage,  threshold) 

[m,  n]  -size  (  inputimage) 

inputimage-real ( inputimage) ; 

for  1-1 :m; 
for  3-l:n; 

if  real{inputiraage(i,j))>255; 

input ima ge{ i, j )-255 ; 

elseif  real(inputimage{i, }))<0; 

inputirnage{  i,  J  )-0; 

end 

end 

end 

F-fftshift(fft2(fftshift( inputimage) ) ) ; 

load  log_oos 

AO-F, *G0; 

A1-F.*G1; 

A2-F,*G2; 

A3-F. *G3 ; 

A4-F, *G4 ; 

A5-F.*G5; 

al-(fftshift(ifft2(fftshift{Al)))); 
a2-(fftshift(ifft2(fftshift(A2)))); 
a3-( f ft8hift( if ft2 { fftshift(A3 ) ) ) ) ; 
a4-(fftshift(ifft2(fftshlft(A4)))); 
a5- ( f ftshift( if ft2{fftshift(A5) ) ) ) ; 

Ll-AO; 

Ll»( fftshift{ if ft2{ fftshift(Ll) ) ) ) ; 

L2-A1+L1; 

12- (fftshift{ifft2(fftshift(L2)))); 

L3-A2+L2; 

13- {fftshlft{ifft2(fftshift(L3)))); 

L4-A3+L3; 

14- (f ftshift ( if ft2(f ftshift{L4) ) ) ) ; 

L5-A4+L4; 

15- (fftshift(ifft2(fftshift(L5) ))); 

cl-(al) ,/ll; 
c2-(a2) ,/12; 
c3-(a3)./13; 
c4-(a4) ./14; 
c5-(a5) ./15; 

for  u-1 :m; 
for  V " 1 : n ; 

if  abs(cl(u, V) ) <-threshoid( 1) ; 
finall ( u, V) -0 ; 
else ; 

finall(u, v)-(al{u, V) ) ; 
end  ; 

if  abs (c2 (u, V) )<-threshold(2) ; 

final2(u,v)-0; 

else; 

final2  (u,  v)-(a2  (u,  V) ) ; 
end; 
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if  abs (c3 (u, V) )<-threshold( 3 ) ; 

finals  (u,  v)-0; 

else; 

final]  (u,v)-(a3(u,v)); 
end ; 

if  abs (c4 (u, V) )<-threshold( 4 ) ; 

final4{u,v)-0; 

else; 

final4  (u  ,  V]  -{a4  (u  ,  v] )  ; 
end; 

if  abs{c5{u,v) ]<-thre9hold(5) ; 

final5(u,v)-0; 

else; 

final5(u,v)-(a5{u,v) )  ; 
end; 

end 

end 


output  image-  ( ll+finall+final2-!-final3+final4+final5 ) ; 


Appendix  B.  Enhancement  Algorithms 


This  appendix  contains  the  code  used  for  the  various  not-wavelet  enhancement  algo¬ 
rithms  outlined  in  this  thesis. 


function  newpicture-adap_filter( picture) 

picture -pad 2 (picture , zeros (21)); 

(A, B] -size (picture) ; 
lovpass -zeros ( A, B ) ; 

inverse-l/(  21''2 ) ; 

for  x-11: (A-11) 
for  y-11: (B-11) 

lowpas3(x, y) -9ura( 9um(picture(x-10 ix+l0,y-10:y+10) ) ) 'Inverse ; 

end 

end 

highpass-5*{ picture -lowpass) ; 

for  x-1 ; A- 20 
for  y-l;B-20 

lowpass2(x, y)-lowpass(x+10,y+10); 
highpass2{x, y)-hlghpass(x+10, y+10) ; 
end 
end 

newpicture-highpass2+lowpass2; 


iThis  function  will  apply  the  pre*emphasis  filter  to  the  Inputlmage. 

function  output  image-emphasis ( inputimage) 

[m,n)-si2e(inputimage} ; 

input  image-real ( inputimage) ; 

for  i-l;m; 
for  3-1 ; n; 

if  real( inputimage( i, 3 ) )>255 ; 

input image( i, j )-255; 

elseif  real(inputimage(i,3))<0: 

inputimage( 1, j )-0; 

end 

end 

end 

F-f ftshif t ( ff t2 ( fftshift( inputimage) ) ) ; 

load  log_cos 

A0-F.*G0; 

A1-F.»G1; 

A2-F.*G2; 

A3-F.*G3; 

A4-F.*G4; 

A5-F.*G5; 

a4-5* ( f ftshif t( ifft2(fftshift(A4 ) ) ) )  ; 

Ll-AO; 

11- ( f ftshift{ ifft2( fftshift(Ll) ) ) ) ; 

L2-A1+L1; 

12- (fftshift(ifft2(fftshift(L2)))); 

L3-A2+L2; 

13-  (fftshift (if ft2(f ftshif t (L3))) ) ; 

L4-A3+L3; 

14- (fftshift(ifft2(fftshift(L4)))); 
outputimage-(a4+14) ; 

output  iinage-heq(  output  image,  255  ) ; 
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Appendix  C.  Matlab  Implementation  of  Biortho  gonal  Wavelets 

C.l  Coefficients 

To  decompose  an  image,  a  scaling  function  and  a  wavelet  function  are  needed.  With 
the  biorthogonal  wavelets  used  in  this  thesis,  the  scaling  function  and  wavelet  function  for 
decomposition  were  defined  in  the  file  FBI.m  as: 


scaletilde=[0  .037828455506995  -.023849465019380... 
-.11062440441842  .37740285561265  .85269867900940... 
.37740285561265  -.11062440441842  -.023849465019380... 
.037828455506995  0] ; 

wavelettilde=[0  0  0  .065438882628938  -.040689417609558... 
-.41809227322221  .78858561640566  -.41809227322221... 
-.040689417609558  .065438882628938  0]  ; 


The  waveforms  generated  and  their  relationship  are  shown  in  figure  27.  The  extra  zeros 
in  the  definitations  are  for  Matlab  purposes,  since  convolution  was  used.  In  addition  to  the 
code  shown  above,  the  following  was  also  in  the  file  FBI.m: 


L-length(scaletilde) ; 
for  n-2 ; L'l; 

wavelet (n-1)- ( ( ) *scaletilde(L-n+l) ; 
end 

wavelet (L-l)-(  ( - 1 ) '‘n)  *scaletllde(L- n) ; 
wavelet(L)-(  ( -  l)''n)  *scaletilde(L-n)  ; 

L-length(wavelet) ; 

for  n-l:L*l; 

scale{n+l)-(('l)^(n'l)) *wavelettilde{L'n+l) ; 
end 

scale  ( 1 )  -  ( ( - 1 )  ■'n)  *wavelettilde(L-n) ; 


When  executed  in  Matlab,  the  values  in  the  vectors  wavelet  and  scale  were  then: 

wavelet  = 

Columns  1  through  7 
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0.0378 


0.0238  -0.1106  -0.3774 


0.8527 


-0.3774  -0.1106 


Columns  8  through  11 

0.0238  0.0378  0 


0 


scale  = 

Columns  1  through  7 

0  0  -0.0654  -0.0407  0.4181  0.7886 

Columns  8  through  11 


-0.0407  -0.0654  0  0 


0.4181 


These  vectors  and  their  relationship  are  shown  in  figure  3 1 .  Once  again,  the  added  zeros 
are  due  to  Matlab  implementation. 


C.2  Decomposition 

With  suitable  coefficients,  decomposition  was  a  matter  of  careful  implementation  of  the 
procedure  outlined  by  Mallat  (18)  and  Laing  (14). 

In  the  file  Down2.m  was  the  function  Down2(inputface,  level).  The  entire  function 
code  was: 


%This  ra-file  will  decompose  a  two  dimensional  image  according  to  % 
IMallat's  method  on  page  685  of  his  1979  wavelet  paper,  % 
Icorrected  by  Laing,  with  inputs  from  Maj.  Greg  Warhola.  I 

Variables 


%  % 

\  % 

I  FACE . image  to  be  decomposed  % 

%  face . dummy  variable  so  FACE  will  not  be  affected  % 

%  LEVEL . how  many  LEVELS  user  wishes  to  decompose  % 

I  set . which  wavelet  set  to  use 

%  taps . how  many  Daubechies  taps  to  use  % 

%  scaletilde . low  frequency  Impulse  response  used  for  decomposition  % 

I  wavelettilde . high  frequency  impulse  response  used  for  decomposition  % 

I  L,K . size  of  image  at  any  given  decomposition  level  % 

I  decompose_level . current  level  of  decomposition  % 

%  gface . face  convolved  with  the  column  vector  wavelettilde'  % 

%  hface . face  convolved  with  the  column  vector  scaletilde'  % 

%  gfacedown . downsampled  gface  (every  other  row,  size  L/2  x  K)  I 

%  hfacedown . downsampled  hface  (every  other  row,  size  L/2  x  K)  % 
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%  ggface . gfacedown  convolved  with  row  vector  wavelettllde  \ 

\  hgface . hfacedown  convolved  with  row  vector  wavelettllde  I 

I  ghface . gfacedown  convolved  with  row  vector  scaletllde  % 

\  hhface . hfacedown  convolved  with  row  vector  scaletilde  % 

%  D3_  . downsampled  ggface  (every  other  column,  size  L/2  x  K/2)  % 

%  D2_  . downsampled  ghface  (every  other  column,  size  L/2  x  K/2)  I 

I  Dl_  . downsampled  hgface  (every  other  column,  size  L/2  x  K/2)  I 

\  A_  . r. . .  .downsampled  hhface  (every  other  column,  size  L/2  x  K/2)  % 

I  LEVEL. . name  of  file  containing  D3_,  D2_,  Dl_,  and  A_P  I 

%  % 


function  Down2 (FACE, LEVEL) ; 

%  This  line  would  let  a  specific  wavelet  set  be  used.  For  this  research, 

%  the  FBI  set  was  used  for  enhancement. 

« 

%x«input( 'Wavelet  set  to  use:  l)Mallat  2)FBI  3)8attle  4)Daubechies. . . ' ) ; 

x-2; 

y-20; 

save  which  x 

if  X— 1 

Mallat; 

M“length(acaletilde) ; 
middle-ceil(M/2) t 
elaelf  x— 2 
FBI; 

M-length ( scaletilde ) ; 
mlddle-cell(M/2); 
elaeif  X--3 
Battle; 

M-length(scaletilde} ; 
mlddle»ceil(M/2) ; 
elaelf  x— 4 

\  y«input( 'How  many  taps?  Enter:  4  6  8  10  12  14  16  18  20...'); 
save  whether  y 
Daubs ; 

H-length(scaletilde) ; 
mlddle-ceil(H/2) ; 
end 

% 

%  The  mirrorface  function  causes  FACE  to  be  quadrupled  to  minimize  edge  effects.  The 
%  original  image  matrix  can  be  thought  of  as  quadrant  one,  quadrant  two  is  the  mirror 
t  image  of  quadrant  one  while  quadrants  three  and  four  are  the  mirror  image  of  quadrants 
%  one  and  two.  This  causes  the  image  to  be  2*L  x  2*K  in  size. 

t 

face*mirrorfaoe(PACE) ; 

Ifaee-FACE; 

[L,Kl-8ize(face) ; 

M-length(scaletilde) ; 

middle-ceil(M/2) ; 

decompose.level-l; 

while  decompoae_level<-LEVEL 

deconq>ose_level  ; 

clear  hface  gface 
for  k-l:K; 
for  1-1:L; 
hface(l,k)-0; 
gface(l,)c)-0; 
for  m-l:M; 

% 

%  arg  is  a  virtual  argument  for  circular  convolution.  The  following 
\  insures  circular  convolution  occurs  regardless  of  the  size  of  the 
I  current  matrix 

t 

arg-l-m+middle; 

while  (arg>L) 1 (arg<-0) 

if  arg>L 

arg-arg-L; 

elaeif  arg<-0 

arg-L+arg; 

end 

end 

I 

%  The  following  is  the  actual  filtering  of  the  scale  and  wavelet 
I  functions  in  the  direction  of  decomposition 

I 
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hface(l,  k)-hface{l,  )c)+f ace(arg,  k)*scaletilde(ni)  ; 

gface(l, k)-gface(l, k)+face(arg, k)*wavelettilde(m); 

end 

end 

end 

clear  gfacedown  hfacedown  ggface  hgface  ghface  hhface 

% 

I  Downsampling  the  columns 

I 

for  l-l:L/2; 

gfacedown(l, : }-gface(2*l, : ) ; 
hfacedown(l, : )-hface(2*l, : )  ; 
end 

for  l-l:L/2; 

for  k-l:K; 

ggface(L/2,K)-0; 

hgface(L/2,K)-0; 

ghface(L/2,K)-0; 

hhface(L/2,K)-0; 

for 

% 

%  arg  is  a  virtual  argument  for  circular  convolution.  The  following 
%  insures  circular  convolution  occurs  regardless  of  the  size  of  the 
\  current  matrix 

t 

arg-k-m+middle; 

while  (arg>K) I (arg<-0) 

if  arg>K 

arg-arg-K; 

elseif  arg<*0 

arg-K+arg; 

end 

end 

% 

1  The  following  is  the  actual  Altering  with  the  scale  and  wavelet 
t  functions  in  the  direction  of  decomposition 

% 

hhface(l,k)-hhface(l,k}+hfacedown(l,arg)*scaletllde(m) ; 

ghface(l,k}-ghface{l,k)+gfacedown(l,arg)*8calettlde(m); 

hgface(l,k)-hgface(l,k)+hfacedown(l,arg)*wavelettllde(m); 

ggface{l,k)-ggface{l,k)+gfaccdown(l,arg)*wavelettllde(m); 

end 

end 

end 

I  Downsampling  the  rows 

\ 

for  k-l:K/2j 

eval(('D3_'  lnt28tr(decompo8e_level)  ' ( wh)-ggface( ; ,2*k) ; ' 1 ) 
eval(('D2_'  lnt28tr(decompo8e_level)  '  ( : ,k)-ghface( : , 2*k) ; ' J ) 
eval(('Dl_'  lnt28tr(decompo8e_level)  ' ( j ,k)-hgface{ : , 2*k); ' ] ) 
eval({'A_'  int28tr(decomp08e_level)  '  { : ,  k)-hhface( : , 2*k) ; ' ] ) 
end 

I 

%  Saving  the  current  decomposition  level  coefficients 

% 

eval{ ( 'save  LEh^EL_’  int2str(decompo8e_level)  '  D3_'  lnt2str(decompo8e_level) 
'  D2_'  int28tr{decompo8e_level)'  Dl_'  int2str(decompo8e_level)  '  A_'  ... 
int2str(decompose_level) ] ) 

eval( I ' face“A_'  int2str(decompo8e_level) 

l-l/2; 

K-K/2; 

decompoae_level-decompo8e_level+l; 

end 
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C.3  Reconstruction 


The  reconstruction  code  was  in  the  file  Up2.m  and  was 


kThis  m-file  will  reconstruct  a  two  dimensional  image  according  to  % 
iMallat's  method  on  page  685  of  his  1979  wavelet  paper,  % 
tcorrected  by  Laing,  with  inputs  from  Maj.  Greg  Warhola.  % 

Variables 


I 

I 

LEVEL . how  many  LEVELS  were  decomposed  % 

set . which  wavelet  set  was  used  % 

taps . how  many  Daubechies  taps  were  used  t 

scale . low  frequency  impulse  response  used  for  decomposition  I 

wavelet . high  frequency  impulse  response  used  for  decomposition  % 

A . current  image  matrix  to  be  upsampled  and  convolved  t 

L,K . size  of  image  at  any  given  reconstruction  level  I 

reconstruct_level . current  level  of  decomposition  % 

LE?VEL_ . name  of  file  containing  D3_,  D2_,  Dl^,  and  A_  % 

D3_ . diagonal  detail  coefficients  at  current  level  \ 

D2_ . horizontal  detail  coefficients  at  current  level  % 

Dl_ . vertical  detail  coefficients  at  current  level  % 

A_ . approximation  coefficients  at  current  level  % 

ggfacetemp.  .  .  . . .D3_  with  rows  upsampled  and  zero  padded  (size  L  x  2*k)  t 

hgfacetemp . D2_  with  rows  upsampled  and  zero  padded  (size  L  x  2*k)  % 

ghfacetemp . Dl_  with  rows  upsampled  and  zero  padded  (size  L  x  2*k)  % 

hhfacetemp . A  with  rows  upsampled  and  zero  padded  (size  L  x  2*k)  % 

ggface . ggfacetemp  convolved  with  the  row  vector  wavelet  % 

hgface . hgfacetemp  convolved  with  the  row  vector  scale  % 

ghface . ghfacetemp  convolved  with  the  row  vector  wavelet  I 

hhface . hhfacetemp  convolved  with  the  row  vector  scale  \ 

gfacedown . ggface+hgface  % 

hfacedown . ghface+hhface  % 

gfacetemp . gfacedown,  cola  upsampled  (size  2*L  x  2*k)  % 

hfacetemp . hfacedown,  cols  upsampled  (size  2*L  x  2*k)  t 

gface . gfacetemp  convolved  with  the  column  vector  wavelet'  I 

hface . hfacetemp  convolved  with  the  column  vector  scale'  % 

face . resulting  reconstructed  image  % 

I 


'unction  face«Gp2 ( LEVEL ) ; 

I  For  this  research,  x-2  was  used  for  the  FBI  wavelet  set. 


load  which 

if  X--1 
Mallat; 
galn*2; 

M“length(8cale) ; 
middle-ceil(M/2) ; 
elseif  X— 2 
FBI; 
gain*l; 

M"length( scale) ; 
middle-cell(M/2) ; 
elseif  X— 3 
Battle; 
gain-1, ■ 

M«length(scale) ; 
mlddle-ceil (M/2 ) ; 
elseif  X— 4 
load  whether 
Daubs; 
gain-1; 

M-length ( scale ) ; 

middle-ceil(M/2); 

end 

reconstruct_level-LEVEL ; 
while  reconstruct_level>0 
reconstruct_level ; 

eval(['load  LEVEL.'  int2str(reconstruct_level)l) 

if  reconstruct.level— LEVEL 
eval(I'A-A_'  int2str(reconstruct_level) 
else 
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A-face; 

end 

M-length(scale)  ; 

[L, K]-size(A) ; 

L-2*L; 

K-2*K; 

middle-ceil(M/2) ; 
enda-floor(M/2 ) ; 

ggfacetemp-zeros(L/2,K) ; 
hgfaceteinp-zeros(L/2,K)  ; 
ghfacetemp-zeros(L/2,K) ; 
hhfaceteinp-zeros(L/2,K) ; 

% 

%  Upsampling  the  rows 

% 

for  k-l!K/2; 

eval( ( 'ggfacetemp( : , 2*k)"D3_'  int2str{recQnatruct_level)  ' ( : *k) ; ' ) ) 
eval( ( 'hgfacetemp( : , 2*k)"D2_'  int2atr(recQnatruct_level)  ' ( s »k) ; ' ] ) 
eval{ [ 'ghfacetemp{ : , 2*k)-Dl_'  int2str{reconatruct_level)  ' ( J »k) ; ' ] ) 
hhfacetemp( : / 2*k)“A( : ,k) ; 
end 

for  l-l!L/2; 

for  k-l:K; 

ggface(l,k)»0; 

hgface(l,k}-0; 

ghface(l,k)*0; 

hhface(l,k)-0; 

for  m-1 sM; 

% 

I  arg  is  a  virtual  argument  for  circular  convolution.  The  following 
I  insures  circular  convolution  occurs  regardless  of  the  size  of  the 
\  current  matrix 

\ 

arg-k-m+mlddle; 

while  (arg>K) I (arg<"0) 

if  arg>K 

arg«arg*K; 

elseif  arg<-0 

arg-K+arg; 

end 

end 

I 

t  The  following  is  the  actual  filtering  of  the  scale  and  wavelet 
I  functions  in  the  direction  of  decomposition 

t 

ggf ace ( 1, k ) -ggf ace( l,k)+ggfacetemp(l, arg) ‘wavelet (ra); 

hgf ace ( 1, k ) -hgf ace ( 1, k) +hgfacetemp(l, arg) ‘scale (m); 

ghf ace ( l,k)-ghf ace (l,k)+ghfacetemp(l, arg) ‘wavelet (m); 

hhface(l/k)>hhface(l,k)+hhfacetemp{l,arg)*scale(m); 

end 

end 

end 

gfacedown-ggface-fhgface; 

hfacedown-ghface+hhface; 

clear  gfacetemp  hfacetemp  gface  hface 

gf acetemp-zeros ( L ; K ) ; 
hfacetemp-zeros(L>K) ; 

t 

k  Upsampling  the  columns 

I 

for  l-l:L/2; 

gfacetemp{2*l, : )-gfacedown{l, : ) ; 
hfacetemp ( 2 ‘1 , : ) -hf acedown ( 1 , : ) ; 
end 

clear  hface  gface 
for  k“l:K; 
for  1-1:L; 
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hface(l,lc)-0; 
gface(l,  )c)-0; 
for 

I 

\  arg  is  a  virtual  argument  for  circular  convolution.  The  following 
\  insures  circular  convolution  occurs  regardless  of  the  size  of  the 
I  current  matrix 

% 

arg-l-m+middle; 

while  (arg>L) I (arg<-0) 

if  arg>L 

arg-arg-Lf 

elseif  arg<-0 

arg-L+arg; 

end 

end 

% 

I  The  following  is  the  actual  filtering  with  the  scale  and  wavelet 
\  functions  in  the  direction  of  decomposition 

\ 

hf  ace  { 1 ,  )c )  -hf  ace  ( 1 .  k )  +gf  acetemp  ( arg ,  k )  ‘wavelet  (m) ; 

gf ace ( 1, k ) -gf ace ( 1, k )+hfacetemp( arg, k) ‘scale (m); 

end 

end 

end 

face- (gface+hface) ‘gain; 

recon8truct_level-recon8truct_level-l; 

end 

I 

k  This  step  ensures  an  image  the  same  size  as  the  original.  (Remember, 

I  the  mirrorimage  function  was  used  to  minimize  edge  effects) 

« 

face-(face(l:L/2,l:K/2)+face(L: •l:L/2+l,l:K/2)+face(l:L/2,K:-l:K/2+l)+face(L:-l:L/2+l,Kt-l:K/2+l))/4; 
return 
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Appendix  D.  Wavelet  Simulation  Code 


The  Matlab  code  in  this  appendix  explains  how  the  simulation  of  visual  perception 
based  on  contrast  sensitivity  using  bi-orthogonal  wavelets  was  implemented. 


\ 

\  This  function  will  create  a  contrast  image  at  each  level  corresponding  to 
I  measured  points  on  the  contrast  threshold  curve,  and  then  use  that  data  to 
t  manipulate  the  detail  coefficients  at  each  level.  The  returned  image,  waveface 
I  will  be  a  wavelet  simulation  of  what  a  person  with  a  given  contrast  threshold 
%  sees  when  looking  at  the  original  image,  FACE, 
t 


% 

% 

VARIABLES 

I 

I  FACE . input  image 

%  face . dummy  variable  so  FACE  will  not  be  changed 

%  threshold,  t . contrast  threshold  used  for  given  person 

I  LEVEL . maximum  level  of  decomposition,  image  slxe  dependant 

I  L,K . size  of  face 

%  LEVEL. . set  of  coefficients  at  current  decomposition  level 

I  D3. . diagonal  details 

%  D2. . horizontal  details 

%  Dl. . vertical  details 

I  A. . approximation  image  at  that  level 

%  C,  dummy . contrast  image  at  that  level 


function  wavef aee-scQtowave( FACE, threshold, seoLEVEL); 

face-FACE; 

t-threshold; 

tL,K]-8ize(face); 
for  1-1:1; 
for  k-l:K; 

if  real(face(l,k))>255; 

faced, k)-255; 

elself  real(face(l, k) }<0; 

face(l,k)-0; 

end 

end 

end 

Down2scowave(face,  t,  scoLETVEL) ; 
waveface-Up2scowave(scoLE:vEL)  ; 


iThis  m-file  will  decompose  a  two  dimensional  image  according  to  % 
iMallat's  method  on  page  685  of  his  1979  wavelet  paper,  % 
^corrected  by  Laing,  with  inputs  from  Kaj .  Greg  Warhola.  I 

Variables 


%  « 

%  I 

\  FACE . image  to  be  decomposed  \ 

\  face . dummy  variable  so  FACE  will  not  be  affected  % 

%  scoLEVEL . how  many  LEVELS  user  wishes  to  decompose  % 

\  set . which  wavelet  set  to  use 

\  taps . . how  many  Daubechies  taps  to  use  % 

t  scaletilde . low  frequency  in^ulse  response  used  for  decomposition  I 

%  scowavelettilde . high  frequency  Impulse  response  used  for  decomposition  % 

\  L,K . size  of  image  at  any  given  decomposition  level  % 

%  decompose.level . current  level  of  decomposition  % 

%  gface . face  convolved  with  the  column  vector  scowavelettilde'  % 

%  hface . face  convolved  with  the  column  vector  scaletilde'  % 

I  gfacedown . downsampled  gface  (every  other  row,  size  L/2  x  K)  \ 

I  hfacedown . downsampled  hface  (every  other  row,  size  L/2  x  K)  % 

%  ggface . gfacedown  convolved  with  row  vector  scowavelettilde  % 

1  hgface . hfacedown  convolved  with  row  vector  scowavelettilde  \ 

%  ghface . gfacedown  convolved  with  row  vector  scaletilde  % 
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I  hhface . hfacedown  convolved  with  row  vector  scaletilde  % 

I  D3^ . downsampled  ggface  (every  other  column,  size  L/2  x  K/2)  I 

\  D2_ . downsampled  ghface  (every  other  column,  size  L/2  x  K/2)  I 

%  Dl_  . downsampled  hgface  (every  other  column,  size  L/2  x  K/2)  I 

I  A_ . downsampled  hhface  (every  other  column,  size  L/2  x  K/2)  % 

\  scoLEr/EL_ . name  of  file  containing  acoD3_,  scoD2_,  scoDl_,  and  scoA_ 

\  t 


function  Down2scowave(FACE,  threshold,  scoLBVEL); 

%  This  line  would  let  a  specific  wavelet  set  be  used.  For  this  research, 

I  the  FBI  set  was  used  for  simulation. 

% 

%x-input( 'Wavelet  set  to  use:  l)Mallat  2)FBI  3)Battle  4)Daubechle3 . . . ' )  ; 
x-3; 

%y-20; 

save  which  x 

if  x— 1 
Mallat; 

M“length(scaletilde} ; 
middle“ceil(M/2) ; 
elseif  X--2 
FBI; 

M-length(scaletilde) ; 
middle-ceil(M/2) j 
elseif  X“*3 
Battle; 

M-length(scaletilde) ; 
middle-ceil(M/2); 
elseif  X— 4 

1  y*input( 'How  many  taps?  Enter:  4  6  0  10  12  14  16  18  20...'); 

save  whether  y 

Daubs; 

M-length(3caletllde) ; 
middle*ceil(H/2) ; 
end 

\  The  nirrorface  function  causes  FACE  to  be  quadrupled  to  minimize  edge  effects.  The 
I  original  image  matrix  can  be  thought  of  as  quadrant  one,  quadrant  two  is  the  mirror 
I  image  of  quadrant  one  while  quadrants  three  and  four  are  the  mirror  Image  of  quadrants 
\  one  and  two.  This  causes  the  image  to  be  2*L  x  2*K  in  size. 

\ 

tL,Kl-size(FACE); 

face-mlrrorface(PACE) ; 

(L,K]*8ize(face) ; 

M-length(scaletllde) ; 

middle-ceil(M/2) ; 

decompoae_levcl-l ; 

while  decompose_level<-8CoLBVEL 

clear  hface  gface  ggface  hgface  ghface  hhface 

i-scoLEVEL*decompo8e_level-l; 

cpl-(2l)+(2*(i'l)); 

if  cpi  >  32 
for  Jc-1:K; 
for  1-1 :L; 
hface(l,)c)-0; 
for  m-l:M; 

t 

\  arg  is  a  virtual  argument  for  circular  convolution.  The  following 
I  insures  circular  convolution  occurs  regardless  of  the  size  of  the 
I  current  matrix 

t 

arg-l-m+middle; 

while  (arg>L) I (arg<-0) 

if  arg>L 

arg-arg-L; 

elseif  arg<-0 

arg-L+arg; 

end 

end 
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% 

%  The  following  is  the  actual  filtering  of  the  scale  and  wavelet 
I  functions  in  the  direction  of  decomposition 

« 

hface{l, k)-hface( 1, k)+face{arg, k) *8caletilde(m) ; 

end 

end 

end 

else 

for  k"l:K; 
for  1-1:L; 
hface(l,k)-0; 
gface(l,k}-0; 
for  ra-1 :M; 

i 

t  arg  is  a  virtual  argument  for  circular  convolution.  The  following 
%  insures  circular  convolution  occurs  regardless  of  the  size  of  the 
%  current  matrix 

\ 

arg-l-m+middlej 

while  (arg>L) 1 (arg<«0) 

if  arg>L 

arg-arg-L; 

elseif  arg<>0 

arg-L+arg; 

end 

end 

4 

4  The  following  is  the  actual  filtering  of  the  scale  and  wavelet 
4  functions  in  the  direction  of  decomposition 

4 

4444444444444444444444444444444444444444444444444444444444444444444 

hfaced,  k)-hface(l,  k)+face(arg,  k)*scaletllde(m) ; 

gface(l,k)«gface(l,k)'<>face(arg/k)*8cowavelettilde(m); 

end 

end 

end 

end 

if  cpi  >  32 
for  1-1 :L; 
for  k-liK; 
hhface{l,k)-0; 

for  m-l:M; 

4444444444444444444444444444444444444444444444444444444444444444444 

4 

4  arg  Is  a  virtual  argument  for  circular  convolution.  The  following 
4  insures  circular  convolution  occurs  regardless  of  the  size  of  the 
4  current  matrix 

4 

4444444444444444444444444444444444444444444444444444444444444444444 

arg-k*m+middle; 

while  (arg>K) i {arg<-0) 

if  arg>K 

arg-arg-K; 

elseif  arg<-0 

arg-K+arg; 

end 

end 

4444444444444444444444444444444444444444444444444444444444444444444 

4 

4  The  following  is  the  actual  filtering  with  the  scale 
4  functions  in  the  direction  of  decomposition 

4 

4444444444444444444444444444444444444444444444444444444444444444444 

hhface(l,  k)-hhface(l,  k)-fhface(l,arg)-*scaletilde(m) ; 

end 

end 

end 

else 

for  1-1:L; 

for  k-l:K; 

ggface(l,k)-0; 

hgface(l,k)-0; 

ghface(l,k)-0; 

hhface(l,k)-0; 

for  m-1 :M; 

4444444444444444444444444444444444444444444444444444444444444444444 

4 

4  arg  is  a  virtual  argument  for  circular  convolution.  The  following 
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%  insures  circular  convolution  occurs  regardless  of  the  size  of  the 
I  current  matrix 

t 

arg-k-m+middle; 

while  (arg>K) I (arg<-0) 

if  arg>K 

arg-arg*K; 

elseif  arg<"0 

arg-K+arg; 

end 

end 

\ 

I  The  following  is  the  actual  filtering  with  the  scale  and  wavelet 
%  functions  in  the  direction  of  decomposition 

% 

ghface(l,  )c)-ghface(l, k)+gface(l, arg)*8caletllde(m) ; 

hgface(l,k)-hgface(l,k)+hface(l,arg)*scowavelettilde(m); 

ggface(l,k)-ggface(l,k)+gface(l,arg)*scowavelettilde(m); 

hhface(l, k)*hhface(l, k)+hface(l,arg)*scaletilde{m) ; 

end 

end 

end 

end 

if  cpi  <  32 

Cgh-abs {ghf ace ,/(hhface+. 000000001) ) ? 

Chg-abs (hgf ace. 000000001) ) ; 

Cgg-abs (ggf ace. /(hhface+. 000000001) ); 

[LC,KC]-8ize(Cgh); 
for  1-1 :LC; 
for  k-l:KC; 

lf(Cgh(l,k}  <•  threahold(cpl) ) 
ghfaced,  k)-0; 

end 

if(Chg(l,k)  <-  threshold(cpl) ) 

hgface(l,k)-0; 

end 

if(Cgg(l,k)  <•  threshold(cpl) ) 

ggface(l,k)-0; 

end 

end 

end 

end 

t 

\  Downsampling  the  rows  and  columns 

t 

if  cpi  <  32 
for  l-liL/2; 
for  k-l:K/2; 

if  (2*1+1<-L)  &  (2*k+l<-K) 

eval( ( ' 8CoD3_'  int28tr(decompo8e_level)  '  (1,  k)-ggface(2*l+l, 2*k+l) ; ' ] ) 
eval{ ( ' scoD2_'  int2str(decompo8e_level)  ' (1, k)-ghfacc(2*l+l, 2*k); '  1 ) 
eval{ ( ' 8coDl_'  int28tr(decompo8e_level)  ' {1,  k)-hgface(2*l/ 2*k+l); ' J ) 
eval{ ( ' 8coA_'  int28tr(decorapo8e_level)  '  (1,  k)-hhface(2*l, 2*k) ; ' ] ) 
elseif  (2*1+1<-L)  «  (2»k+l>K) 

eval( ( 'scoD3_'  int2str(decompo8e_level)  '  (1,  k)-ggface(2*l+l/ 1) ; ' ] ) 
eval( ( ' scoD2_'  int2str(decompo8e_level)  '  (1,  k)-ghface{2*l+l,2*k) ; ' ) ) 
eval( [ '8coDl_'  int2str(decompo8e_level)  '  (l,k)-hgface( 2*1,1) ; ' ) ) 
eval( ( 'scoA_'  int28tr(decompo8e_level)  ' (l,k)-hhface(2*l, 2*k) ; ' ) ) 
elseif  (2*1+1>L)  i  (2*k+l<-K) 

eval( [ '8coD3_'  int28tr(decompo8e_level)  '  (l,k)-ggface(l, 2*k+l) ; ' J ) 
eval( [ ' 8coD2_'  int28tr(decompo8e_level)  ' (1, k)-ghface(l, 2*k) ; ' ] ) 
eval( ( ' 8coDl_'  int28tr(decomp08e_level)  '  {1,  k)-hgface{2*l, 2*k+l) ; ' J ) 
eval( [ ' 8coA_'  int28tr(decompo8e_level)  ' (1, k)-hhface(2*l, 2*k) ; ' ) ) 
else 

eval{ t '8CoD3_'  int2str{decompo8e_level)  ' (1' J')“99face(l,l) ; ' ) ) 

eval( ( '8coD2_'  int2atr(decompo8e_level)  ' (l,k)-ghface{l, 2*k) ; ' ) ) 

eval( [ ' 8CoDl_'  lnt28tr(dccompoae_level)  ' (l,k)-hgf ace (2*1,1) ; ' J ) 

eval( [ ' scoA_'  int28tr(decompo8e:.level)  ' (l,k)-hhface(2*l,2*k) ; ' J ) 

end 

end 

end 

else 

eval( { ' scoD3_'  int2str(decompose_level)  '-zeros(L/2,K/2);' J) 
eval( [ 'scoD2_'  lnt2str(decompose_level)  '-zeros( V2,K/2) ; ' J ) 
eval{ [ ' 8coD1_'  int28tr(decompo8e_level)  ’ -zero8{L/2,K/2); ' ] ) 
for  l-l:L/2; 
for  k-l:K/2; 

eval( [ 'scoA_'  int28tr(decomp08e_level)  ' (1, k)-hhface(2*l,2*k) ; ' ] ) 

end 

end 

end 
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% 

%  Saving  the  current  decomposition  level  coefficients 


% 

eval( [ ' save  scoLEVEL_'  int2str(decompo8e_level)  '  scoD3_'  lnt2str{decompose_level)  ... 
’  scoD2_'  int2str(decompose_level)  '  scoDl_'  int2str(decompQse_level)  '  9coA_'  ... 
int2str(decompose_level) 1 ) 

eval( [ ' face-scoA^'  int2str(decompose_level) 

I.-L/2; 

K-K/2 ; 

decompose_level-decompose_level+l; 

end 


^This  m’file  will  reconstruct  a  two  dimensional  image  according  to  % 

%Mallat's  method  on  page  685  of  his  1979  wavelet  paper,  I 
%corrected  by  Laing,  with  inputs  from  Maj .  Greg  Warhola.  % 

Variables 

t  I 

t  t 

%  scoLEJVEL . how  many  LEVELS  were  decomposed  % 

\  set . which  wavelet  set  was  used  % 

I  taps . how  many  Oaubechies  taps  were  used  % 

\  scale . low  frequency  Impulse  response  used  for  decomposition  % 

%  wavelet . high  frequency  impulse  response  used  for  decomposition  I 

\  A . current  image  matrix  to  be  upsampled  and  convolved  % 

\  L,K . size  of  image  at  any  given  reconstruction  level  % 

%  recon8truct_l6vel . current  level  of  decomposition  % 

%  9coLEVEL_ . name  of  file  containing  D3_,  D2_,  Dl.,  and  A_  % 

%  D3_ . diagonal  detail  coefficients  at  current  level  I 

%  D2_ . horizontal  detail  coefficients  at  current  level  % 

I  Dl_ . vertical  detail  coefficients  at  current  level  % 

I  A_ . approximation  coefficients  at  current  level  % 

t  ggfacetemp . 03_  with  rows  upsampled  and  zero  padded  (size  L  x  2*k)  % 

I  hgfacetemp . 02.  with  rows  upsampled  and  zero  padded  (size  L  x  2*k)  % 

\  ghfacetemp . 01.  with  rows  upsampled  and  zero  padded  (size  L  x  2*k)  % 

I  hhfacetemp . A  with  rows  upsampled  and  zero  padded  (size  L  x  2*k)  % 

\  ggface . ggfacetemp  convolved  with  the  row  vector  wavelet  % 

%  hgface . hgfacetemp  convolved  with  the  row  vector  scale  % 

%  ghface . ghfacetemp  convolved  with  the  row  vector  wavelet  % 

%  hhfaee . hhfacetemp  convolved  with  the  row  vector  scale  % 

\  gfacedown . ggface+hgface  \ 

%  hfacedown . ghface-t-hhface  % 

%  gfacetemp . gfacedown,  cols  upsampled  (size  2*L  x  2*k)  % 

\  hfacetemp . hfacedown,  cols  upsampled  (size  2*L  x  2*k)  % 

\  gface . gfacetemp  convolved  with  the  column  vector  wavelet'  \ 

%  hface . hfacetemp  convolved  with  the  column  vector  scale'  \ 

\  face . resulting  reconstructed  image  \ 

I  % 

function  face-Up28cowave(8coLEVEL) ; 

\  For  this  research,  x-2  was  used  for  the  FBI  wavelet  set. 

load  which 

if  x— 1 
Mallat; 
galn-2; 

M-length{8cale) ; 

middle-ceil(M/2); 

elseif  X— 2 

FBI; 

gain-l; 

M-length(8cale); 
mlddlc-ceil(M/2) ; 
elseif  X--3 
Battle; 
gain-1; 

M-length(scale)  ; 
middle-ceil(M/2) ; 
elseif  X— 4 
load  whether 
Daubs; 
gain-1; 

M-length(8cale) ; 
middle-ceil (M/2 ) ; 
end 
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reconstruct_level-scoLEVEL ; 
while  reconstruct_level>0 

clear  hface  gface  ggface  hgface  ghface  hhface 

eval( ( ' load  scoLBVEL_'  int2str<reconatruct_level) ] ) 

if  reconstract_level""SCoLEIVEL 

eval( [ 'A-scoA_'  int2str{recon8truct_level} 

else 

A-face; 

end 

M-length(scale) ; 

[L,K3-8i2e(A) ; 

L-2*L; 

K-2*K; 

middle-cell (M/2) ; 
ends-floor(M/2}  ; 

ggfacetemp-2eros(L/2, K) ; 
hgfacetemp-zeros(L/2i K) ; 
ghfacetemp-zero8(L/2, K) ; 
hhfacetemp-2eroa{L/2, K) ; 

% 

%  Upsampling  the  rows 

% 

slze(ggfaceteinp)  ; 

eval(I '8i2e(8coD3_'  lnt2str(reconatruct_level)  ')»'!) 
for  1c-1:K/2; 

eval( [ 'ggfacetemp( : , 2*k)-8CoD3„'  int28tr(recon8truct_level)  ' ( : ,k); ' J ) 
eval( [ 'hgfacetemp( : , 2»k)-8C0D2_'  int2atr(recon8truct_ievel)  ' ( : .k) ; ' 1 ) 
eval( [ 'ghfacetemp( : , 2*k)«8CoDl_'  lnt28tr(recon8truct_level)  ' ( : / k) ; ' ] ) 
hhfacetemp( s , 2*k)-A( i ,  k) ; 
end 

for  l-l:L/2; 
for  k-l:K; 
ggface(l, k)"0; 
hgface(l,k)-0; 
ghface(l,  k)<*0,' 
hhface(l,k)-0; 

for  m-l:M; 

% 

\  arg  is  a  virtual  argument  for  circular  convolution.  The  following 
t  insures  circular  convolution  occurs  regardless  of  the  size  of  the 
I  current  matrix 

% 

arg-k-m+raiddle; 

while  (arg>K) I (arg<-0} 

if  arg>K 

arg-arg-K; 

elseif  arg<-0 

arg-K+arg; 

end 

end 

% 

\  The  following  is  the  actual  filtering  of  the  scale  and  wavelet 
%  functions  in  the  direction  of  decomposition 

\ 

ggface(l,k}-ggface(l,k)-fggfacetemp(l,arg)*wavelet(m); 

hgface(l,k}-hgface(l,k)+hgfacetemp(l/arg)*8cale(m); 

ghface(l,k)-ghface(l,k)+ghfacetemp(l,arg)*wavelet(m); 

hhface(l,k)-hhface(l,k)+hhfacetemp(l,arg)*8cale(m); 

end 

end 

end 

gfacedown-ggface+hgface; 

hfacedown-ghface+hhface; 

clear  gfacetemp  hfacetemp  gface  hface 

gfacetemp-2ero8(L,K) ; 
hfacetemp-2eroB(L, K) ; 
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I 

%  Upsampling  the  columns 

t 

for  l-l:L/2; 

gfacetemp(2*l, : }-gfacedown(l, : ) ; 
hfacetemp(2*l, : )*hfacedown(l, : ) ; 
end 

clear  hface  gface 
for  k-l:K; 
for  1-1:L; 
hface(l,lc)-0; 
gface{l,k)-0; 
for  m“l:M; 

t 

%  arg  is  a  virtual  argument  for  circular  convolution.  The  following 
%  insures  circular  convolution  occurs  regardless  of  the  size  of  the 
%  current  matrix 

% 

arg-l-m+middle; 

while  (arg>L) I (arg<-0) 

if  arg>L 

arg-arg-L; 

elseif  arg<-0 

arg-L+arg; 

end 

end 

« 

\  The  following  is  the  actual  filtering  with  the  scale  and  wavelet 
I  functions  in  the  direction  of  decomposition 

« 

hface(l,)c)»hface(l,)c)+gfacetemp(arg,k)*wavelet(m); 

gf  ace  ( 1,  Jc ) -gf  ace  (l,)c)+hfacetemp  (arg,  k)  *80810  (m); 

end 

end 

end 

face"(gface‘)‘hface}  *gain; 
reoon8truct_level-reconstruct_level • 1 ; 
end 

\ 

%  This  step  ensures  an  image  the  same  size  as  the  original.  (Remember, 

%  the  mirrorimage  function  was  used  to  minimize  edge  effects) 

\ 

face-(face(l:L/2,l!K/2)+face(L:-l:L/2+l,l:K/2)+face(l:L/2,K:-l:K/2+l)+face(L:‘ia/2+l,K:*l:K/2+l)}/4; 

return 
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